Something to get started to create an anisotropic conductivity tensor for S4L from solving a Poisson equation with boundary conditions (running an EM simulation and using its potential). Once the field is created in the analysis tab it can be used as a part of a new EM simulation for any material selected to have inhomogeneous and anisotropic properties.
import numpy as np
import s4l_v1.document as document
import s4l_v1.model as model
import s4l_v1.simulation.emlf as emlf
from s4l_v1 import ReleaseVersion, Unit
from s4l_v1.model import Vec3, Translation, Rotation
import s4l_v1.units as units
import s4l_v1.analysis as analysis
import XCoreModeling
def CreateMockAnisotropicConductivity(sim):
# Transversal and longitudinal sigmas
sigma1 = 3.0
sigma2 = 1.0
# here conductivity has 6 component (anisotropic tensor)
results = sim.Results()
p = results["Overall Field"].Outputs["EM Potential(x,y,z,f0)"]
p.Update()
e = results["Overall Field"].Outputs["EM E(x,y,z,f0)"]
e.Update()
e_field = e.Data.Field(0)
nc = e.Data.Grid.NumberOfCells
print "E Field shape", e_field.shape
print "Number of Cells", nc
# E Field shape (nc, 3) - Cell centered
# Calculate sigma on cell centers
print "Calculate anisotropic sigma values"
sigma_wm = np.zeros([nc,6]) # 6 unique values
# tensor component indexing for solver
iu3 = (np.array([0, 1, 2, 0, 1, 0], dtype=np.int64), np.array([0, 1, 2, 1, 2, 2], dtype=np.int64))
for idx, val in enumerate(e_field):
tmp = np.zeros([3,3], dtype=np.complex)
tmp[0] = val
if (not np.isnan(val).any() and np.inner(val, val) != 0):
P = np.dot(tmp.T, tmp) / np.inner(val,val)
else:
P = np.full((3, 3), 0, dtype=np.complex)
sigma_tensor = sigma1*P + sigma2*( np.identity(P.shape[0]) - P )
sigma_wm[idx] = sigma_tensor[iu3].real
sigma_swapped = sigma_wm[:,[0,1,2,3,5,4]] # Fix bug where yz and xz were swapped
tensor_conductivity = analysis.core.DoubleFieldData()
tensor_conductivity.NumberOfSnapshots = 1
tensor_conductivity.NumberOfComponents = 6
tensor_conductivity.Grid = e.Data.Grid.Clone()
tensor_conductivity.ValueLocation = analysis.core.eValueLocation.kCellCenter
tensor_conductivity.SetField(0, sigma_swapped) #sigma_wm) # Fix bug where yz and xz were swapped
tensor_conductivity.Quantity.Name = 'Anisotropic Conductivity From E Field'
tensor_conductivity.Quantity.Unit = analysis.core.Unit('S/m')
producer_t = analysis.core.TrivialProducer()
producer_t.SetDataObject(tensor_conductivity)
producer_t.Name = 'Anisotropic Conductivity'
document.AllAlgorithms.Add(producer_t)