refactors makeLayer in example

Author: StephenNneji

ratapi/examples/domains/domains_XY_model.py

@@ -1,8 +1,9 @@
 """Custom model file for the domains custom XY example."""
 
-import math
+from math import sqrt
 
 import numpy as np
+from scipy.special import erf
 
 
 def domains_XY_model(params, bulk_in, bulk_out, contrast, domain):
@@ -19,13 +20,13 @@ def domains_XY_model(params, bulk_in, bulk_out, contrast, domain):
     z = np.arange(0, 141)
 
     # Make the volume fraction distribution for our Silicon substrate
-    [vfSilicon, siSurf] = makeLayer(z, -25, 50, 1, subRough, subRough)
+    [vfSilicon, siSurf] = make_layer(z, -25, 50, 1, subRough, subRough)
 
     # ... and the Oxide ...
-    [vfOxide, oxSurface] = makeLayer(z, siSurf, oxideThick, 1, subRough, subRough)
+    [vfOxide, oxSurface] = make_layer(z, siSurf, oxideThick, 1, subRough, subRough)
 
     # ... and also our layer.
-    [vfLayer, laySurface] = makeLayer(z, oxSurface, layerThick, 1, subRough, layerRough)
+    [vfLayer, laySurface] = make_layer(z, oxSurface, layerThick, 1, subRough, layerRough)
 
     # Everything that is not already occupied will be filled will water
     totalVF = vfSilicon + vfOxide + vfLayer
@@ -53,7 +54,7 @@ def domains_XY_model(params, bulk_in, bulk_out, contrast, domain):
     return SLD, subRough
 
 
-def makeLayer(z, prevLaySurf, thickness, height, Sigma_L, Sigma_R):
+def make_layer(z, prevLaySurf, thickness, height, Sigma_L, Sigma_R):
     """Produce a layer, with a defined thickness, height and roughness.
 
     Each side of the layer has its own roughness value.
@@ -63,12 +64,9 @@ def makeLayer(z, prevLaySurf, thickness, height, Sigma_L, Sigma_R):
     right = prevLaySurf + thickness
 
     # Make our heaviside
-    a = (z - left) / ((2**0.5) * Sigma_L)
-    b = (z - right) / ((2**0.5) * Sigma_R)
+    erf_left = erf((z - left) / (sqrt(2) * Sigma_L))
+    erf_right = erf((z - right) / (sqrt(2) * Sigma_R))
 
-    erf_a = np.array([math.erf(value) for value in a])
-    erf_b = np.array([math.erf(value) for value in b])
-
-    VF = np.array((height / 2) * (erf_a - erf_b))
+    VF = np.array((0.5 * height) * (erf_left - erf_right))
 
     return VF, right

ratapi/examples/normal_reflectivity/custom_XY_DSPC.py

@@ -52,10 +52,10 @@ def custom_XY_DSPC(params, bulk_in, bulk_out, contrast):
     z = np.arange(0, 141)
 
     # Make our Silicon substrate
-    vfSilicon, siSurf = layer(z, -25, 50, 1, subRough, subRough)
+    vfSilicon, siSurf = make_layer(z, -25, 50, 1, subRough, subRough)
 
     # Add the Oxide
-    vfOxide, oxSurface = layer(z, siSurf, oxideThick, 1, subRough, subRough)
+    vfOxide, oxSurface = make_layer(z, siSurf, oxideThick, 1, subRough, subRough)
 
     # We fill in the water at the end, but there may be a hydration layer between the bilayer and the oxide,
     # so we start the bilayer stack an appropriate distance away
@@ -66,15 +66,15 @@ def custom_XY_DSPC(params, bulk_in, bulk_out, contrast):
     headThick = vHead / lipidAPM
 
     # ... and make a box for the volume fraction (1 for now, we correct for coverage later)
-    vfHeadL, headLSurface = layer(z, watSurface, headThick, 1, bilayerRough, bilayerRough)
+    vfHeadL, headLSurface = make_layer(z, watSurface, headThick, 1, bilayerRough, bilayerRough)
 
     # ... also do the same for the tails
     # We'll make both together, so the thickness will be twice the volume
     tailsThick = (2 * vTail) / lipidAPM
-    vfTails, tailsSurf = layer(z, headLSurface, tailsThick, 1, bilayerRough, bilayerRough)
+    vfTails, tailsSurf = make_layer(z, headLSurface, tailsThick, 1, bilayerRough, bilayerRough)
 
     # Finally the upper head ...
-    vfHeadR, headSurface = layer(z, tailsSurf, headThick, 1, bilayerRough, bilayerRough)
+    vfHeadR, headSurface = make_layer(z, tailsSurf, headThick, 1, bilayerRough, bilayerRough)
 
     # Making the model
     # We've created the volume fraction profiles corresponding to each of the groups.
@@ -120,7 +120,7 @@ def custom_XY_DSPC(params, bulk_in, bulk_out, contrast):
     return SLD, subRough
 
 
-def layer(z, prevLaySurf, thickness, height, Sigma_L, Sigma_R):
+def make_layer(z, prevLaySurf, thickness, height, Sigma_L, Sigma_R):
     """Produce a layer, with a defined thickness, height and roughness.
 
     Each side of the layer has its own roughness value.