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.