Cycles: Fix three numerical issues in the fresnel, normal map and Beckmann code

- In fresnel_dielectric, the differentials calculation sometimes divided by zero.
- When the normal map was (0.5, 0.5, 0.5), the code would try to normalize a zero vector. Now, it just uses the regular normal as a fallback.
- The approximate error function used in Beckmann sampling sometimes overflowed to inf while calculating r^16. The final value is 1 - 1/r^16, however,
  so now it just returns 1 if the computation would overflow otherwise.

intern/cycles/kernel/closure/bsdf_util.h

@@ -80,7 +80,7 @@ ccl_device float fresnel_dielectric(
 		return 1; // total internal reflection
 	}
 	else {
-		float dnp = sqrtf(arg);
+		float dnp = max(sqrtf(arg), 1e-7f);
 		float nK = (neta * cos)- dnp;
 		*T = -(neta * I)+(nK * Nn);
 #ifdef __RAY_DIFFERENTIALS__

intern/cycles/kernel/svm/svm_tex_coord.h

@@ -312,7 +312,7 @@ ccl_device void svm_node_normal_map(KernelGlobals *kg, ShaderData *sd, float *st
 
 		/* apply normal map */
 		float3 B = sign * cross(normal, tangent);
-		N = normalize(color.x * tangent + color.y * B + color.z * normal);
+		N = safe_normalize(color.x * tangent + color.y * B + color.z * normal);
 
 		/* transform to world space */
 		object_normal_transform(kg, sd, &N);
@@ -330,14 +330,18 @@ ccl_device void svm_node_normal_map(KernelGlobals *kg, ShaderData *sd, float *st
 		if(space == NODE_NORMAL_MAP_OBJECT || space == NODE_NORMAL_MAP_BLENDER_OBJECT)
 			object_normal_transform(kg, sd, &N);
 		else
-			N = normalize(N);
+			N = safe_normalize(N);
 	}
 
 	float strength = stack_load_float(stack, strength_offset);
 
 	if(strength != 1.0f) {
 		strength = max(strength, 0.0f);
-		N = normalize(ccl_fetch(sd, N) + (N - ccl_fetch(sd, N))*strength);
+		N = safe_normalize(ccl_fetch(sd, N) + (N - ccl_fetch(sd, N))*strength);
+	}
+
+	if(N == make_float3(0.0f, 0.0f, 0.0f)) {
+		N = ccl_fetch(sd, N);
 	}
 
 	stack_store_float3(stack, normal_offset, N);

intern/cycles/util/util_math_fast.h

@@ -547,6 +547,9 @@ ccl_device_inline float fast_erff(float x)
 	const float a5 = 0.0002765672f;
 	const float a6 = 0.0000430638f;
 	const float a = fabsf(x);
+	if(a >= 12.3f) {
+		return copysignf(1.0f, x);
+	}
 	const float b = 1.0f - (1.0f - a);  /* Crush denormals. */
 	const float r = madd(madd(madd(madd(madd(madd(a6, b, a5), b, a4), b, a3), b, a2), b, a1), b, 1.0f);
 	const float s = r * r;  /* ^2 */