State and Materials

The State class plays a central role in every simulation. It’s main purpose is to keep track of the current simulation state, such as the magentization and the time. Moreover, the state object holds the mesh along with domain information and helps to connect the different solver modules through virtual attributes.

Domain Management

In the simplest case, a state object is generated with a plain mesh object

mesh  = Mesh("sample.xml")
state = State(mesh)

Typical structures simulated with magnum.fe have a size of nanometers to microns. Since magnum.fe expects all quantities to be given in SI units, the mesh size should be given in meters. However, meshes with dimensions of \(1 \cdot 10^{-9}\) lead to problems both at mesh generation and solution of finite-element systems. Hence it is recommended to use nanometers as mesh scale and set a proper scaling in the State initializer.

mesh  = Mesh("sample.xml")
state = State(mesh, scale = 1e-9)

There are two possibilties to define domains in the mesh. Either the domain information is embedded directly in the mesh XML file or the domain markers are provided in separate MeshFunction objects. For details, see Meshes and domains.

Domains in the mesh files are addressed by integer IDs. The State class can be initialized with a domain name mapping for convenient access A more convenient way to address domains is by

In order to address domains by given names the State class can be provided with a mapping of domain names to IDs. This significantly increases the readability of simulation scripts. Moreover, some domain names are of special importance since different magnum.fe module solve problems on certain regions, e.g. the LLG is by default solved on the region named ‘magnetic’. Naming domains is done on intialization of the state object.

state = State(mesh,
  cell_domains = {
    'magnetic':   (2,4),
    'conducting': (1,2,3,4,5),
    'fixed'     : 2,
    'free'      : 4
  },
  facet_domains = {
    'contact1': 1,
    'contact2': 2
  },
  scale = 1e-9)

The naming of cell domains (volumes) and facet domains (areas) is defined by distinct Python dictionaries. The given names can be used in a variety of methods of the State class.

# integrate function f over 'fixed' domain
f = Function(state.FunctionSpace())
assemble(f * state.dx('fixed'))

# Compute volume of 'free' domain
V = state.volume('free')

# Compute area of 'contact1' domain
A = state.area('contact1')

Besides these generic uses of the domain names almost every solver class (e.g. field terms, time integrators) can be configured to solve on certain domains. For instance the LLGCvode class solves the LLG on the ‘magnetic’ region be default, but can easily be configured to solve on another region. The same applies to the effective field contributions.

# solve the LLG with exchange field only on 'some_region'
exchange = ExchangeField(region = 'some_region')
llg = LLGCvode([exchange], region = 'some_region')

Furthermore, nonlocal field terms such as the demagnetization field and the Oersted field can be configured with a source and a target region.

# solve the demag field generated by region 'a' on region 'b'
demag = DemagField(source_region = 'a', target_region = 'b')

This feature is for instance useful for the simulation of multilayer structures with a free and fixed magnetic layer that acts as spin polarizer. By restricting the source region to the free layer, the demagnetization generated by the fixed layer can be neglected while still considering its magnetic properties for the computation of the spin accumulation.

When using the FEM/BEM method for the solution of the nonlocal field terms, see Open Boundary Problems, an additional compute_region can be defined. This region has to contain both the source and the target region and defines the area that is solved with pure FEM. In the case of multilayer structures with alternating magnetic and nonmagnetic layers, setting the compute region to ‘all’ usually reduces the boundary that has to be considered for the boundary element method. Furthermore, by an appropriate choice of the compute region, small distances of boundary elements, that lead to inaccurate BEM coefficients, can be avoided. For further information refer to the reference of the DemagField and OerstedField class.

Material Definition

Material parameters are managed by the State class and made accessible via the material attribute. They can be set either globally or in a specific domain.

# global assignment of multiple parameters
state.material = Material(
  alpha = 1.0,
  ms    = 8e5
)

# global assignment of single paramter
state.material.alpha = 0.1

# assign multiple parameters in region named 'fixed'
state.material['fixed'] = Material(
  K_uni      = 1e6,
  K_uni_axis = (1, 0, 0)
)

# assign single parameter in region named 'free'
state.material['free'].K_uni = 0.0

# assign spatially varying parameter
state.material['magnetic'].K_uni = Expression("x[0]")

Depending on the applied modules, a multitude of material parameters has to be set. Internally, all material parameters are discretized by elementwise constant functions. The parameters are applied in the order of appearance in the script, so the subsequent call of

state.material.alpha = 1.0
state.material.alpha = 0.1

will lead to alpha = 0.1. The attribute getters of state.material return a FEniCS Function object, so proper assignment can be checked by writing the material function to a file.

File("alpha.pvd") << state.material.alpha

Attributes

The state object is used to communicate any simulation state between the different modules and classes of magnum.fe. Different modules require the state to provide different attributes. For instance the DemagField module requires the state object to provide the magnetization state.m, while the OerstedField module depends on the electric current state.j.

Setting an attribute on the state object is as easy as calling state.j = Constant((1.0, 0.0, 0.0)). Note that attributes can be either simple scalars or GenericFunction objects, see FEniCS manual. GenericFunction is a superclass to Constant, Expression and Function objects. Any GenericFunction object will be automatically converted to a Function.

# Each of these commands will result in state.m being a Function object
state.m = Constant((1, 0, 0))
state.m = Expression(("x[0]", "0", "0"))
state.m = Function(...)

GenericFunction object are not only converted to Function objects, but also decorated with some useful methods that are aware of the named domains.

# Return function m restricted on the 'fixed' region
state.m['fixed']

# Compute componentwise average of m
state.m.average()

# Compute componentwise average of m on 'free' region
state.m.average('free')

Virtual Attributes

Once set, simple attributes as introduced in the preceding section remain constant. However, attributes are often required to change depending on the simulation state. For instance, in case of an AC current, the electric current obviously depends on the time. In order to provide attributes depending on other attributes, magnum.fe provides the concept of virtual attributes. Virtual attribues can be defined by a function that takes the state as single parameter and returns either a scalar or a GenericFunction.

# Define an AC current with amplitude A and frequency freq
state.j = lambda state: Constant((A * sin(2*pi*freq*state.t), 0, 0))

The main purpose of virtual attributes is to provide certain attributes to certain solvers, e.g. the OerstedField requires state.j to be defined. However, you may define arbitrary virtual attributes for your own convenience.

While virtual attributes are very useful to define simple time dependencies, e.g. for the definition of an AC current as shown above, they might also be used to describe complicated dependencies of the attribute from the complete simulation state. For instance the electric current might depend on a potential difference between two contacts and additionally on the magnetization configuration when dealing with spin diffusion models. magnum.fe offers a number of predefined functions that can be used as virtual attributes. Simply speaking virtual attributes provide an interface for the flexible combination of different models.

# Compute currect j with spin diffusion model
spin_acc = SpinAccumulationForPotential()
state.j = spin_acc.j()

# Compute Oersted field
oersted = OerstedField()
oersted.h(state)

# Alternatively register Oersted field computation as virtual attribute
state.h = oersted.h()

Sometimes multiple modules will require the same virtual attribute in a single time step. If the computation of the virtual attribute is expensive, the virtual attribute should only be computed once and the result should be reused as long as possible. For this purpose virtual attribues provide a simple caching mechanism. If the virtual attributes depends on one are more other attributes, caching can be enabled by returning not the result itself, but a tuple containing the results as first entry and the names of the depending attributes as further entries.

# compute energy -inner(m, h) and recompute only on change of state.m
state.E = lambda state: (assemble(-inner(state.m, h)*dx), "m")

# set m
state.m = Constant((1.0, 0.0, 0.0))
energy = state.E # trigger computation of E
energy = state.E # takes E from cache

# reset m
state.m = Constant((0.0, 1.0, 0.0))
energy = state.E # trigger computation of E

Predefined virtual attributes of magnum.fe usually implement reasonable caching so you dont’ have to worry about that.

Masked State

Complex micromagnetic models might consist of a multitude of regions and some subproblems, i.e. partial differential equations, might be defined only on certains regions. For this use case, magnum.fe provides the possibility to mask a state to a certain domain. A masked state acts exactly like a regular state, except that the mesh and all function valued attributes are region to the masked region.

# Create a masked state for region 'magnetic'
masked_state = state['magnetic']

# Override values of attribute j only in 'magnetic' region
masked_state.j = Constant((1.0, 0.0, 0.0))

# This returns True
assemble(Constant(1.)*state.dx('magnetic')) == assemble(Contant(1.)*masked_state.dx())

Class Reference

State

class magnumfe.State(mesh, cell_domains={}, facet_domains={}, celldomains={}, facetdomains={}, cell_markers=None, facet_markers=None, material=None, maxwell_material=None, scale=1.0, t=0.0, check_mesh=True, **kwargs)

This class holds the complete state of the simulation and provides some convenience wrappers for the handling of named domains. Furthermore it provides an interface for attribute handling and caching.

Domain/Material Examples
state = State(mesh,
  cell_domains  = {'magnetic': (1,3), 'conducting': 2, 'iron': 1, 'cobalt': 3},
  facet_domains = {'contact1': 1, 'contact2': 2},
  m = Expression(...)
)

# Set materials for different regions
state.material['iron']       = Material(...)
state.material['cobalt']     = Material(...)
state.material['conducting'] = Material(...)

# Use integration measures with named domains
assemble(Constant(1.0) * state.dx('all'))       # All named domains
assemble(Constant(1.0) * state.dx('magnetic'))  # Magnetic region
assemble(Constant(1.0) * state.dx('!magnetic')) # Nonmagnetic region
assemble(Constant(1.0) * state.dx(1))           # Region by ID

# Compute average of magnetization m
state.m.average()       # Over whole space
state.m.average('iron') # Over iron region

# Crop magnetization to subdomain and save as PVD
File("m_iron.pvd") << state.m['iron']

# Normalize the magnetization
state.m.normalize()
Attribute Examples
# initialize state with constant magnetization. m is automatically
# interpolated on a suitable discrete function space.
state = State(mesh,
  m = Constant((1.0, 0.0, 0.0)),
)

# define a current as function of the time
state.j = lambda state: Constant((state.t * 1e12, 0.0, 0.0))

# define some functional depending on m
state.E = lambda state: (assemble(inner(state.m, state.m)*state.dx()), "m")

# dependencies can be arbirtrarily nested
state.E_times_2 = lambda state: (2*state.E, "E")

# all values are cached according to their dependencies
state.E_times_2  # triggers computation of E_times_2
state.E_times_2  # taken from cache, no computation

state.m = Constant((0.0, 1.0, 0.0))

state.E_times_2  # triggers computation of E_times_2
Arguments
mesh (dolfin.Mesh)
The mesh including all subdomains as dolfin.MeshDomains
cell_domains (dict)
naming of the cell subdomains, at least the subdomains ‘magnetic’ and ‘conducting’ should be defined.
facet_domains (dict)
naming of the facet subdomains
cell_markers (dolfin.MeshFunction)
cell domain markers as MeshFunction if domains are not embedded in mesh
facet_markers (dolfin.MeshFunction)
facet domain markers as MeshFunction if domains are not embedded in mesh
material (Material)
the material of the sample. If material differs from subdomain to subdomain, use material setters instead.
scale (float)
the spatial scaling of the mesh. Use 1e-9 if you use nanometers as length measure.
t (float)
the time
check_mesh (bool)
Check mesh quality and print warning if quality is bad
**kwargs (dict)
add any state variables like magnetization (m) or spin diffusion (s). Expressions are automatically interpolated on the corresponding discrete spaces.
FunctionAssigner(family='Lagrange', degree=1, rank=1)

Returns a function assigner

Arguments
family (str)
element type of the function space: “CG”, “DG”, ...
degree (int)
the degree of the function space (default=1)
rank (int)
the rank of the target function space (0 for scalar space, 1 for vector space(default))
Example
assigner_V_to_VV = state.FunctionAssigner()
assigner_VV_to_V = state.FunctionAssigner(rank=0)
Returns
dolfin.FunctionAssigner
the function assigner
FunctionSpace(family='Lagrange', degree=1, rank=0)

Returns a function space for a mesh.

Arguments
family (str)
element type of the function space
degree (int)
the degree of the function space (default=1)
rank (int)
the rank of the function space (0 for scalar space, 1 for vector space)
Example
V = state.FunctionSpace("CG", 2)
VV = state.FunctionSpace("CG", 1, rank=1)
Returns
dolfin.FunctionSpace
the function space
M_inv_diag(domain='all', weight=None)

Returns the inverse lumped mass matrix for the vector function space. The result is cached.

Arguments
domain (str)
the domain on which the mass matrix is restricted
weight (str)
weight to be used for assembly, either name of material attribute (e.g. “ms”) or None
Returns
dolfin.Matrix
the matrix
VectorFunctionSpace(family='Lagrange', degree=1)

Returns the vector-function space for a mesh.

Arguments
family (str)
element type of the function space
degree (int)
the degree of the function space (default=1)
Returns
dolfin.FunctionSpace
the function space
area(domain='all')

Computes the area of a certain facet domain.

Arguments
domain (string / int)
domain identifier (either string or id)
Returns
float
the area of the region
dP(domain='all', domaintype='cell', graceful=False)

Convenience wrapper for integral-point measure. If the mesh does not contain any cell domains, the measure for the whole mesh is returned.

Arguments
domain (string / int)
name or ID of domain
domaintype (string)
either ‘cell’ or ‘facet’
Returns
dolfin:Measure
the measure
dS(domain='all', dS=None, graceful=False)

Convenience wrapper for integral-interior-facet measure. If the mesh does not contain any cell domains, the measure for the whole mesh is returned.

Arguments
domain (string / int)
name or ID of domain
dS (dolfin:Measure)
alternative measure
Returns
dolfin:Measure
the measure
domain_ids(domain, domaintype='cell')

Get the domain IDs for a named domain.

Arguments
domain (string/ int)
name or ID of domain
domaintype (string)
domain type, (cell, facet)
Returns
[int]
List of domain IDs
ds(domain='all', ds=None, intersect=None, graceful=False)

Convenience wrapper for integral-facet measure. If the mesh does not contain any cell domains, the measure for the whole mesh is returned.

Arguments
domain (string / int)
name or ID of domain
ds (dolfin:Measure)
alternative measure
intersect (string / int)
name or ID of anothe domain to intersect with
Returns
dolfin:Measure
the measure
dx(domain='all', dx=None, graceful=False)

Convenience wrapper for integral-cell measure. If the mesh does not contain any cell domains, the measure for the whole mesh is returned.

Arguments
domain (string / int)
name or ID of domain
dx (dolfin:Measure)
alternative measure
Returns
dolfin:Measure

the measure

# Compute volume of magnetic domain assemble(Constant(1.0) * state.dx(‘magnetic’))

# Compute volume of domain with ID 3 assemble(Constant(1.0) * state.dx(3))

get_or_create(key, func)

Simple cache helper that checks for a cache hit in the generic state.cache dict and populates the cache with the result of func evaluation otherwise.

Arguments
key (str)
cache key
func (str)
function call to be cached
interpolate(expressions, family='Lagrange', degree=1)

Interpolates a collection of expressions, each defined for a specific domain, on a suitable function space. Also works for a single expression.

Arguments
expressions (dict)
Expressions to be interpolated.
Returns
dolfin.Function
Interpolated function
Example
f1 = state.interpolate({
  'magnetic':  Constant((1.0, 0.0, 0.0)),
  '!magnetic': Constant((0.0, 1.0, 0.0))
})

f2 = state.interpolate(Constant(1.0))
intersect(*domains, **kwargs)

Creates the intersection of multiple domains and returns the IDs.

Examples
# compute the intersection volume of the regions 'fixed' and 'free'
state.volume(state.intersect('magnetic', 'free'))

# compute the intersection area of the regions 'contact1' and 'contact2'
state.area(state.intersect('contact1', 'contact2'))
Arguments
domains (string/ int)
name or ID of domain
domaintype (string)
type of domain (‘cell’ or ‘facet’, default is ‘cell’)
Returns
the intersection region
mask(domain)

Return masked state on a certain region.

Argument
domain
str
domain name
Returns
MaskedState
the masked state
step(integrators, dt)

Calls the step method on every integrator and increases t.

Arguments
integrators (list or Integrator)
Integrator on which the step method is called. Can be either a list of integrators or a single integrator.
union(*domains, **kwargs)

Creates the union of multiple domains and returns the IDs.

Examples
# compute the total volume of the regions 'fixed' and 'free'
state.volume(state.union('fixed', 'free'))

# compute the total area of the regions 'contact1' and 'contact2'
state.area(state.union('contact1', 'contact2', domaintype='facet'))
Arguments
domains (string/ int)
name or ID of domain
domaintype (string)
type of domain (‘cell’ or ‘facet’, default is ‘cell’)
Returns
the union region
update(*names)

Visits the given attributes in order to update them if needed. This function should be called if references to attributes are cached.

Arguments
*names (string)
list of attribute names
uuid(*names)

Returns a uuid for a given list of attributes. The uuid changes if the attribute if reset.

Arguments
*names (string)
list of attribute names
Returns
the uuid
volume(domain='all')

Computes the volume of a certain cell domain.

Arguments
domain (string / int)
domain identifier (either string or id)
Returns
float
the volume of the region

Material

class magnumfe.Material(*args, **kwargs)

Create a material.

Example
# create simple material
a = Material(alpha = 1.0, ms = 8e5)

# create material with material a as template
b = Material(a, alpha = 0.1)
Arguments
*args
templates
**kwargs
material constants
parameters()

Returns all material parameters as dict.

Returns
dict
the parameters
static py()

Return an instance of Material with the material properties similar to Permalloy.