API Reference¶

CASHOCS is a computational, adjoint based shape optimization and optimal control software for python.

CASHOCS is based on the finite element package FEniCS and uses its high-level unified form language UFL to treat general PDE constrained optimization problems, in particular, shape optimization and optimal control problems.

Here, we detail the (public) API of cashocs.

For a more hands-on approach, we recommend the tutorial, which shows many examples from PDE constrained optimization that can be efficiently solved with CASHOCS.

PDE Constrained Optimization Problems¶

If you are using CASHOCS to solve PDE constrained optimization problems, you should use the following two classes, for either optimal control or shape optimization problems.

OptimalControlProblem¶

class cashocs.OptimalControlProblem(state_forms, bcs_list, cost_functional_form, states, controls, adjoints, config=None, riesz_scalar_products=None, control_constraints=None, initial_guess=None, ksp_options=None, adjoint_ksp_options=None, scalar_tracking_forms=None, min_max_terms=None, desired_weights=None)¶

Bases: cashocs._interfaces.optimization_problem.OptimizationProblem

Implements an optimal control problem.

This class is used to define an optimal control problem, and also to solve it subsequently. For a detailed documentation, see the examples in the tutorial. For easier input, when considering single (state or control) variables, these do not have to be wrapped into a list. Note, that in the case of multiple variables these have to be grouped into ordered lists, where state_forms, bcs_list, states, adjoints have to have the same order (i.e. [y1, y2] and [p1, p2], where p1 is the adjoint of y1 and so on.

This is used to generate all classes and functionalities. First ensures consistent input, afterwards, the solution algorithm is initialized.

Parameters

state_forms (ufl.form.Form or list[ufl.form.Form]) – The weak form of the state equation (user implemented). Can be either a single UFL form, or a (ordered) list of UFL forms.
bcs_list (list[fenics.DirichletBC] or list[list[fenics.DirichletBC]] or fenics.DirichletBC or None) – The list of DirichletBC objects describing Dirichlet (essential) boundary conditions. If this is None, then no Dirichlet boundary conditions are imposed.
cost_functional_form (ufl.form.Form or list[ufl.form.Form]) – UFL form of the cost functional.
states (fenics.Function or list[fenics.Function]) – The state variable(s), can either be a fenics.Function, or a list of these.
controls (fenics.Function or list[fenics.Function]) – The control variable(s), can either be a fenics.Function, or a list of these.
adjoints (fenics.Function or list[fenics.Function]) – The adjoint variable(s), can either be a fenics.Function, or a (ordered) list of these.
config (configparser.ConfigParser or None, optional) – The config file for the problem, generated via cashocs.create_config(). Alternatively, this can also be None, in which case the default configurations are used, except for the optimization algorithm. This has then to be specified in the solve method. The default is None.
riesz_scalar_products (None or ufl.form.Form or list[ufl.form.Form], optional) – The scalar products of the control space. Can either be None, a single UFL form, or a (ordered) list of UFL forms. If None, the \(L^2(\Omega)\) product is used. (default is None).
control_constraints (None or list[fenics.Function] or list[float] or list[list[fenics.Function]] or list[list[float]], optional) – Box constraints posed on the control, None means that there are none (default is None). The (inner) lists should contain two elements of the form [u_a, u_b], where u_a is the lower, and u_b the upper bound.
initial_guess (list[fenics.Function], optional) – List of functions that act as initial guess for the state variables, should be valid input for fenics.assign(). Defaults to None, which means a zero initial guess.
ksp_options (list[list[str]] or list[list[list[str]]] or None, optional) – A list of strings corresponding to command line options for PETSc, used to solve the state systems. If this is None, then the direct solver mumps is used (default is None).
adjoint_ksp_options (list[list[str]] or list[list[list[str]]] or None, optional) – A list of strings corresponding to command line options for PETSc, used to solve the adjoint systems. If this is None, then the same options as for the state systems are used (default is None).
scalar_tracking_forms (dict or list[dict] or None, optional) – A list of dictionaries that define scalar tracking type cost functionals, where an integral value should be brought to a desired value. Each dict needs to have the keys 'integrand' and 'tracking_goal'. Default is None, i.e., no scalar tracking terms are considered.
min_max_terms (dict or list[dict] or None, optional) – Additional terms for the cost functional, not to be used directly.
desired_weights (list[float] or None, optional) – A list of values for scaling the cost functional terms. If this is supplied, the cost functional has to be given as list of summands. The individual terms are then scaled, so that term i has the magnitude of desired_weights[i] for the initial iteration. In case that desired_weights is None, no scaling is performed. Default is None.

Return type

OptimalControlProblem

Examples

Examples how to use this class can be found in the tutorial.

compute_adjoint_variables()¶

Solves the adjoint system.

This can be used for debugging purposes and solver validation. Updates / overwrites the user input for the adjoint variables. The solve of the corresponding state system needed to determine the adjoints is carried out automatically.

Returns
Return type: None

compute_gradient()¶

Solves the Riesz problem to determine the gradient.

This can be used for debugging, or code validation. The necessary solutions of the state and adjoint systems are carried out automatically.

Returns: A list consisting of the (components) of the gradient.
Return type: list[fenics.Function]

compute_state_variables()¶

Solves the state system.

This can be used for debugging purposes and to validate the solver. Updates and overwrites the user input for the state variables.

Returns
Return type: None

gradient_test(u=None, h=None, rng=None)¶

Taylor test to verify that the computed gradient is correct for optimal control problems.

Parameters

u (list[fenics.Function] or None, optional) – The point, at which the gradient shall be verified. If this is None, then the current controls of the optimization problem are used. Default is None.
h (list[fenics.Function] or None, optional) – The direction(s) for the directional (Gateaux) derivative. If this is None, one random direction is chosen. Default is None.
rng (numpy.random.RandomState or None, optional) – A numpy random state for calculating a random direction

Returns

The convergence order from the Taylor test. If this is (approximately) 2 or larger,: everything works as expected.

Return type

float

inject_post_hook(function)¶

Changes the a-posteriori hook of the OptimizationProblem

Parameters: function (function) – A custom function without arguments, which will be called after the computation of the gradient(s)
Returns
Return type: None

inject_pre_hook(function)¶

Changes the a-priori hook of the OptimizationProblem

Parameters: function (function) – A custom function without arguments, which will be called before each solve of the state system
Returns
Return type: None

inject_pre_post_hook(pre_function, post_function)¶

Changes the a-priori (pre) and a-posteriori (post) hook of the OptimizationProblem

Parameters

pre_function (function) – A function without arguments, which is to be called before each solve of the state system
post_function (function) – A function without arguments, which is to be called after each computation of the (shape) gradient

Returns

Return type

None

solve(algorithm=None, rtol=None, atol=None, max_iter=None)¶

Solves the optimization problem by the method specified in the config file.

Updates / overwrites states, controls, and adjoints according to the optimization method, i.e., the user-input fenics.Function() objects.

Parameters

algorithm (str or None, optional) – Selects the optimization algorithm. Valid choices are 'gradient_descent' or 'gd' for a gradient descent method, 'conjugate_gradient', 'nonlinear_cg', 'ncg' or 'cg' for nonlinear conjugate gradient methods, 'lbfgs' or 'bfgs' for limited memory BFGS methods, 'newton' for a truncated Newton method, and 'pdas' or 'primal_dual_active_set' for a primal dual active set method. This overwrites the value specified in the config file. If this is None, then the value in the config file is used. Default is None.
rtol (float or None, optional) – The relative tolerance used for the termination criterion. Overwrites the value specified in the config file. If this is None, the value from the config file is taken. Default is None.
atol (float or None, optional) – The absolute tolerance used for the termination criterion. Overwrites the value specified in the config file. If this is None, the value from the config file is taken. Default is None.
max_iter (int or None, optional) – The maximum number of iterations the optimization algorithm can carry out before it is terminated. Overwrites the value specified in the config file. If this is None, the value from the config file is taken. Default is None.

Returns

Return type

None

Notes

If either rtol or atol are specified as arguments to the solve call, the termination criterion changes to:

a purely relative one (if only rtol is specified), i.e.,

\[|| \nabla J(u_k) || \leq \texttt{rtol} || \nabla J(u_0) ||.\]

a purely absolute one (if only atol is specified), i.e.,

\[|| \nabla J(u_k) || \leq \texttt{atol}.\]

a combined one if both rtol and atol are specified, i.e.,

\[|| \nabla J(u_k) || \leq \texttt{atol} + \texttt{rtol} || \nabla J(u_0) ||.\]

supply_adjoint_forms(adjoint_forms, adjoint_bcs_list)¶

Overwrites the computed weak forms of the adjoint system.

This allows the user to specify their own weak forms of the problems and to use cashocs merely as a solver for solving the optimization problems.

Parameters

adjoint_forms (ufl.Form or list[ufl.Form]) – The UFL forms of the adjoint system(s).
adjoint_bcs_list (list[fenics.DirichletBC] or list[list[fenics.DirichletBC]] or fenics.DirichletBC or None) – The list of Dirichlet boundary conditions for the adjoint system(s).

Returns

Return type

None

supply_custom_forms(derivatives, adjoint_forms, adjoint_bcs_list)¶

Overrides both adjoint system and derivatives with user input.

This allows the user to specify both the derivatives of the reduced cost functional and the corresponding adjoint system, and allows them to use cashocs as a solver.

ShapeOptimizationProblem¶

class cashocs.ShapeOptimizationProblem(state_forms, bcs_list, cost_functional_form, states, adjoints, boundaries, config=None, shape_scalar_product=None, initial_guess=None, ksp_options=None, adjoint_ksp_options=None, scalar_tracking_forms=None, min_max_terms=None, desired_weights=None)¶

Bases: cashocs._interfaces.optimization_problem.OptimizationProblem

A shape optimization problem.

This class is used to define a shape optimization problem, and to solve it subsequently. For a detailed documentation, we refer to the tutorial. For easier input, when consider single (state or control) variables, these do not have to be wrapped into a list. Note, that in the case of multiple variables these have to be grouped into ordered lists, where state_forms, bcs_list, states, adjoints have to have the same order (i.e. [y1, y2] and [p1, p2], where p1 is the adjoint of y1 and so on.

Parameters

state_forms (ufl.Form or list[ufl.Form]) – The weak form of the state equation (user implemented). Can be either a single UFL form, or a (ordered) list of UFL forms.
bcs_list (list[fenics.DirichletBC] or list[list[fenics.DirichletBC]] or fenics.DirichletBC or None) – The list of DirichletBC objects describing Dirichlet (essential) boundary conditions. If this is None, then no Dirichlet boundary conditions are imposed.
cost_functional_form (ufl.Form or list[ufl.Form]) – UFL form of the cost functional.
states (fenics.Function or list[fenics.Function]) – The state variable(s), can either be a fenics.Function, or a list of these.
adjoints (fenics.Function or list[fenics.Function]) – The adjoint variable(s), can either be a fenics.Function, or a (ordered) list of these.
boundaries (fenics.MeshFunction) – fenics.MeshFunction that indicates the boundary markers.
config (configparser.ConfigParser or None, optional) – The config file for the problem, generated via cashocs.create_config(). Alternatively, this can also be None, in which case the default configurations are used, except for the optimization algorithm. This has then to be specified in the solve method. The default is None.
shape_scalar_product (ufl.form.Form or None, optional) – The bilinear form for computing the shape gradient (or gradient deformation). This has to use fenics.TrialFunction and fenics.TestFunction objects to define the weak form, which have to be in a fenics.VectorFunctionSpace of continuous, linear Lagrange finite elements. Moreover, this form is required to be symmetric.
initial_guess (list[fenics.Function] or None, optional) – List of functions that act as initial guess for the state variables, should be valid input for fenics.assign(). Defaults to None, which means a zero initial guess.
ksp_options (list[list[str]] or list[list[list[str]]] or None, optional) – A list of strings corresponding to command line options for PETSc, used to solve the state systems. If this is None, then the direct solver mumps is used (default is None).
adjoint_ksp_options (list[list[str]] or list[list[list[str]]] or None) – A list of strings corresponding to command line options for PETSc, used to solve the adjoint systems. If this is None, then the same options as for the state systems are used (default is None).
scalar_tracking_forms (dict or list[dict] or None, optional) – A list of dictionaries that define scalar tracking type cost functionals, where an integral value should be brought to a desired value. Each dict needs to have the keys 'integrand' and 'tracking_goal'. Default is None, i.e., no scalar tracking terms are considered.
min_max_terms (dict or list[dict] or None, optional) – Additional terms for the cost functional, not to be used directly.
desired_weights (list[float] or None, optional) – A list of values for scaling the cost functional terms. If this is supplied, the cost functional has to be given as list of summands. The individual terms are then scaled, so that term i has the magnitude of desired_weights[i] for the initial iteration. In case that desired_weights is None, no scaling is performed. Default is None.

Return type

ShapeOptimizationProblem

compute_adjoint_variables()¶

Solves the adjoint system.

This can be used for debugging purposes and solver validation. Updates / overwrites the user input for the adjoint variables. The solve of the corresponding state system needed to determine the adjoints is carried out automatically.

Returns
Return type: None

compute_shape_gradient()¶

Solves the Riesz problem to determine the shape gradient.

This can be used for debugging, or code validation. The necessary solutions of the state and adjoint systems are carried out automatically.

Returns: The shape gradient.
Return type: fenics.Function

compute_state_variables()¶

Solves the state system.

This can be used for debugging purposes and to validate the solver. Updates and overwrites the user input for the state variables.

Returns
Return type: None

get_vector_field()¶

Returns the TestFunction for defining shape derivatives.

Additionally constrained problems¶

ConstrainedOptimalControlProblem¶

class cashocs.ConstrainedOptimalControlProblem(state_forms, bcs_list, cost_functional_form, states, controls, adjoints, constraints, config=None, riesz_scalar_products=None, control_constraints=None, initial_guess=None, ksp_options=None, adjoint_ksp_options=None, scalar_tracking_forms=None)¶

Bases: cashocs._constraints.constrained_problems.ConstrainedOptimizationProblem

Parameters

state_forms (ufl.form.Form or list[ufl.form.Form]) – The weak form of the state equation (user implemented). Can be either a single UFL form, or a (ordered) list of UFL forms.
bcs_list (list[fenics.DirichletBC] or list[list[fenics.DirichletBC]] or fenics.DirichletBC or None) – The list of DirichletBC objects describing Dirichlet (essential) boundary conditions. If this is None, then no Dirichlet boundary conditions are imposed.
cost_functional_form (ufl.form.Form or list[ufl.form.Form]) – UFL form of the cost functional.
states (fenics.Function or list[fenics.Function]) – The state variable(s), can either be a fenics.Function, or a list of these.
controls (fenics.Function or list[fenics.Function]) – The control variable(s), can either be a fenics.Function, or a list of these.
adjoints (fenics.Function or list[fenics.Function]) – The adjoint variable(s), can either be a fenics.Function, or a (ordered) list of these.
constraints (EqualityConstraint or InequalityConstraint or list[EqualityConstraint, InequalityConstraint]) – (A list of) additional equality and inequality constraints for the problem.
config (configparser.ConfigParser or None, optional) – The config file for the problem, generated via cashocs.create_config(). Alternatively, this can also be None, in which case the default configurations are used, except for the optimization algorithm. This has then to be specified in the solve method. The default is None.
riesz_scalar_products (None or ufl.form.Form or list[ufl.form.Form], optional) – The scalar products of the control space. Can either be None, a single UFL form, or a (ordered) list of UFL forms. If None, the \(L^2(\Omega)\) product is used. (default is None).
control_constraints (None or list[fenics.Function] or list[float] or list[list[fenics.Function]] or list[list[float]], optional) – Box constraints posed on the control, None means that there are none (default is None). The (inner) lists should contain two elements of the form [u_a, u_b], where u_a is the lower, and u_b the upper bound.
initial_guess (list[fenics.Function], optional) – List of functions that act as initial guess for the state variables, should be valid input for fenics.assign(). Defaults to None, which means a zero initial guess.
ksp_options (list[list[str]] or list[list[list[str]]] or None, optional) – A list of strings corresponding to command line options for PETSc, used to solve the state systems. If this is None, then the direct solver mumps is used (default is None).
adjoint_ksp_options (list[list[str]] or list[list[list[str]]] or None, optional) – A list of strings corresponding to command line options for PETSc, used to solve the adjoint systems. If this is None, then the same options as for the state systems are used (default is None).
scalar_tracking_forms (dict or list[dict] or None, optional) – A list of dictionaries that define scalar tracking type cost functionals, where an integral value should be brought to a desired value. Each dict needs to have the keys 'integrand' and 'tracking_goal'. Default is None, i.e., no scalar tracking terms are considered.

Return type

None

inject_post_hook(function)¶

Changes the a-posteriori hook of the OptimizationProblem

Parameters: function (function) – A custom function without arguments, which will be called after the computation of the gradient(s)
Returns
Return type: None

inject_pre_hook(function)¶

Changes the a-priori hook of the OptimizationProblem

Parameters: function (function) – A custom function without arguments, which will be called before each solve of the state system
Returns
Return type: None

inject_pre_post_hook(pre_function, post_function)¶

Changes the a-priori (pre) and a-posteriori (post) hook of the OptimizationProblem

Parameters

pre_function (function) – A function without arguments, which is to be called before each solve of the state system
post_function (function) – A function without arguments, which is to be called after each computation of the (shape) gradient

Returns

Return type

None

solve(method='Augmented Lagrangian', tol=0.01, max_iter=25, inner_rtol=None, inner_atol=None, constraint_tol=None, mu_0=None, lambda_0=None)¶

Solves the constrained problem

Parameters

method (str, optional) – The solution algorithm, either an augmented Lagrangian method (“Augmented Lagrangian”, “AL”) or quadratic penalty method (“Quadratic Penalty”, “QP”)
tol (float, optional) – An overall tolerance to be used in the algorithm. This will set the relative tolerance for the inner optimization problems to tol. Default is 1e-2.
max_iter (int, optional) – Maximum number of iterations for the outer solver. Default is 25.
inner_rtol (float or None, optional) – Relative tolerance for the inner problem. Default is None, so that inner_rtol = tol is used.
inner_atol (float or None, optional) – Absolute tolerance for the inner problem. Default is None, so that inner_atol = tol/10 is used.
constraint_tol (float or None, optional) – The tolerance for the constraint violation, which is desired. If this is None (default), then this is specified as tol/10.
mu_0 (float or None, optional) – Initial value of the penalty parameter. Default is None, which means that mu_0 = 1 is used.
lambda_0 (list[float] or None, optional) – Initial guess for the Lagrange multipliers. Default is None, which corresponds to a zero guess.

Returns

Return type

None

total_constraint_violation()¶

Compute the total constraint violation.

Returns: The 2-norm of the total constraint violation.
Return type: float

ConstrainedShapeOptimizationProblem¶

class cashocs.ConstrainedShapeOptimizationProblem(state_forms, bcs_list, cost_functional_form, states, adjoints, boundaries, constraints, config=None, shape_scalar_product=None, initial_guess=None, ksp_options=None, adjoint_ksp_options=None, scalar_tracking_forms=None)¶

Bases: cashocs._constraints.constrained_problems.ConstrainedOptimizationProblem

Parameters

state_forms (ufl.Form or list[ufl.Form]) – The weak form of the state equation (user implemented). Can be either a single UFL form, or a (ordered) list of UFL forms.
bcs_list (list[fenics.DirichletBC] or list[list[fenics.DirichletBC]] or fenics.DirichletBC or None) – The list of DirichletBC objects describing Dirichlet (essential) boundary conditions. If this is None, then no Dirichlet boundary conditions are imposed.
cost_functional_form (ufl.Form or list[ufl.Form]) – UFL form of the cost functional.
states (fenics.Function or list[fenics.Function]) – The state variable(s), can either be a fenics.Function, or a list of these.
adjoints (fenics.Function or list[fenics.Function]) – The adjoint variable(s), can either be a fenics.Function, or a (ordered) list of these.
boundaries (fenics.MeshFunction) – fenics.MeshFunction that indicates the boundary markers.
constraints (EqualityConstraint or InequalityConstraint or list[EqualityConstraint, InequalityConstraint]) – (A list of) additional equality and inequality constraints for the problem.
config (configparser.ConfigParser or None, optional) – The config file for the problem, generated via cashocs.create_config(). Alternatively, this can also be None, in which case the default configurations are used, except for the optimization algorithm. This has then to be specified in the solve method. The default is None.
shape_scalar_product (ufl.form.Form or None, optional) – The bilinear form for computing the shape gradient (or gradient deformation). This has to use fenics.TrialFunction and fenics.TestFunction objects to define the weak form, which have to be in a fenics.VectorFunctionSpace of continuous, linear Lagrange finite elements. Moreover, this form is required to be symmetric.
initial_guess (list[fenics.Function] or None, optional) – List of functions that act as initial guess for the state variables, should be valid input for fenics.assign(). Defaults to None, which means a zero initial guess.
ksp_options (list[list[str]] or list[list[list[str]]] or None, optional) – A list of strings corresponding to command line options for PETSc, used to solve the state systems. If this is None, then the direct solver mumps is used (default is None).
adjoint_ksp_options (list[list[str]] or list[list[list[str]]] or None) – A list of strings corresponding to command line options for PETSc, used to solve the adjoint systems. If this is None, then the same options as for the state systems are used (default is None).
scalar_tracking_forms (Optional[Dict]) –

Return type

None

inject_post_hook(function)¶

Changes the a-posteriori hook of the OptimizationProblem

Parameters: function (function) – A custom function without arguments, which will be called after the computation of the gradient(s)
Returns
Return type: None

inject_pre_hook(function)¶

Changes the a-priori hook of the OptimizationProblem

Parameters: function (function) – A custom function without arguments, which will be called before each solve of the state system
Returns
Return type: None

inject_pre_post_hook(pre_function, post_function)¶

Changes the a-priori (pre) and a-posteriori (post) hook of the OptimizationProblem

Parameters

pre_function (function) – A function without arguments, which is to be called before each solve of the state system
post_function (function) – A function without arguments, which is to be called after each computation of the (shape) gradient

Returns

Return type

None

solve(method='Augmented Lagrangian', tol=0.01, max_iter=25, inner_rtol=None, inner_atol=None, constraint_tol=None, mu_0=None, lambda_0=None)¶

Solves the constrained problem

Parameters

method (str, optional) – The solution algorithm, either an augmented Lagrangian method (“Augmented Lagrangian”, “AL”) or quadratic penalty method (“Quadratic Penalty”, “QP”)
tol (float, optional) – An overall tolerance to be used in the algorithm. This will set the relative tolerance for the inner optimization problems to tol. Default is 1e-2.
max_iter (int, optional) – Maximum number of iterations for the outer solver. Default is 25.
inner_rtol (float or None, optional) – Relative tolerance for the inner problem. Default is None, so that inner_rtol = tol is used.
inner_atol (float or None, optional) – Absolute tolerance for the inner problem. Default is None, so that inner_atol = tol/10 is used.
constraint_tol (float or None, optional) – The tolerance for the constraint violation, which is desired. If this is None (default), then this is specified as tol/10.
mu_0 (float or None, optional) – Initial value of the penalty parameter. Default is None, which means that mu_0 = 1 is used.
lambda_0 (list[float] or None, optional) – Initial guess for the Lagrange multipliers. Default is None, which corresponds to a zero guess.

Returns

Return type

None

total_constraint_violation()¶

Compute the total constraint violation.

Returns: The 2-norm of the total constraint violation.
Return type: float

EqualityConstraint¶

class cashocs.EqualityConstraint(variable_function, target, measure=None)¶

Bases: cashocs._constraints.constraints.Constraint

Models an (additional) equality constraint.

Parameters

variable_function (ufl.Form or ufl.core.expr.Expr) – Either a ufl Form (when we have a scalar / integral constraint) or an ufl expression (when we have a pointwise constraint), which models the part that is to be constrained
target (float) – The target (rhs) of the equality constraint.
measure (fenics.Measure or None, optional) – A measure indicating where a pointwise constraint should be satisfied.

Return type

None

constraint_violation()¶

Computes the constraint violation for self

Returns: The computed violation
Return type: float

InequalityConstraint¶

class cashocs.InequalityConstraint(variable_function, lower_bound=None, upper_bound=None, measure=None)¶

Bases: cashocs._constraints.constraints.Constraint

Parameters

variable_function (ufl.Form or ufl.core.expr.Expr) – Either a ufl Form (when we have a scalar / integral constraint) or an ufl expression (when we have a pointwise constraint), which models the part that is to be constrained
lower_bound (float or fenics.Function or None, optional) – The lower bound for the inequality constraint
upper_bound (float or fenics.Function or None, optional) – The upper bound for the inequality constraint
measure (fenics.Measure or None, optional) – A measure indicating where a pointwise constraint should be satisfied.

constraint_violation()¶

Computes the constraint violation of self

Returns: The computed constraint violation.
Return type: float

Command Line Interface¶

For the command line interface of CASHOCS, we have a mesh conversion tool which converts GMSH .msh files to .xdmf ones, which can be read with the import mesh functionality. It’s usage is detailed in the following.

cashocs-convert¶

Convert GMSH to XDMF.

usage: cashocs-convert [-h] infile outfile

Positional Arguments¶

infile: GMSH file to be converted, has to end in .msh
outfile: XDMF file to which the mesh shall be converted, has to end in .xdmf

MeshQuality¶

class cashocs.MeshQuality¶

A class used to compute the quality of a mesh.

This class implements either a skewness quality measure, one based on the maximum angle of the elements, or one based on the radius ratios. All quality measures have values in \([0,1]\), where 1 corresponds to the reference (optimal) element, and 0 corresponds to degenerate elements.

Examples

This class can be directly used, without any instantiation, as shown here

import cashocs

mesh, _, _, _, _, _ = cashocs.regular_mesh(10)

min_skew = cashocs.MeshQuality.min_skewness(mesh)
avg_skew = cashocs.MeshQuality.avg_skewness(mesh)

min_angle = cashocs.MeshQuality.min_maximum_angle(mesh)
avg_angle = cashocs.MeshQuality.avg_maximum_angle(mesh)

min_rad = cashocs.MeshQuality.min_radius_ratios(mesh)
avg_rad = cashocs.MeshQuality.avg_radius_ratios(mesh)

min_cond = cashocs.MeshQuality.min_condition_number(mesh)
avg_cond = cashocs.MeshQuality.avg_condition_number(mesh)

This works analogously for any mesh compatible with FEniCS.

Return type: None

static avg_condition_number(mesh)¶

Computes average mesh quality based on the condition number of the reference mapping.

This quality criterion uses the condition number (in the Frobenius norm) of the (linear) mapping from the elements of the mesh to the reference element. Computes the average of the condition number over all elements.

Parameters: mesh (fenics.Mesh) – The mesh, whose quality shall be computed.
Returns: The average mesh quality based on the condition number.
Return type: float

classmethod avg_maximum_angle(mesh)¶

Computes the average quality measure based on the largest angle.

This measures the relative distance of a triangle’s angles or a tetrahedron’s dihedral angles to the corresponding optimal angle. The optimal angle is defined as the angle an equilateral (and thus equiangular) element has. This is defined as

\[1 - \max\left( \frac{\alpha - \alpha^*}{\pi - \alpha^*} , 0 \right),\]

where \(\alpha\) is the corresponding (dihedral) angle of the element and \(\alpha^*\) is the corresponding (dihedral) angle of the reference element. To compute the quality measure, the average of this expression over all elements and all of their (dihedral) angles is computed.

Parameters: mesh (fenics.Mesh) – The mesh, whose quality shall be computed.
Returns: The average quality, based on the maximum angle measure.
Return type: float

static avg_radius_ratios(mesh)¶

Computes the average radius ratio of the mesh.

This measures the ratio of the element’s inradius to it’s circumradius, normalized by the geometric dimension. This is computed via

\[d \frac{r}{R},\]

where \(d\) is the spatial dimension, \(r\) is the inradius, and \(R\) is the circumradius. To compute the (global) quality measure, the average of this expression over all elements is returned.

Parameters: mesh (fenics.Mesh) – The mesh, whose quality shall be computed.
Returns: The average radius ratio of the mesh.
Return type: float

classmethod avg_skewness(mesh)¶

Computes the average skewness of the mesh.

This measure the relative distance of a triangle’s angles or a tetrahedron’s dihedral angles to the corresponding optimal angle. The optimal angle is defined as the angle an equilateral, and thus equiangular, element has. The skewness lies in \([0,1]\), where 1 corresponds to the case of an optimal (equilateral) element, and 0 corresponds to a degenerate element. The skewness corresponding to some (dihedral) angle \(\alpha\) is defined as

\[1 - \max \left( \frac{\alpha - \alpha^*}{\pi - \alpha*} , \frac{\alpha^* - \alpha}{\alpha^* - 0} \right),\]

where \(\alpha^*\) is the corresponding angle of the reference element. To compute the quality measure, the average of this expression over all elements and all of their (dihedral) angles is computed.

Parameters: mesh (fenics.Mesh) – The mesh, whose quality shall be computed.
Returns: The average skewness of the mesh.
Return type: float

static min_condition_number(mesh)¶

Computes minimal mesh quality based on the condition number of the reference mapping.

This quality criterion uses the condition number (in the Frobenius norm) of the (linear) mapping from the elements of the mesh to the reference element. Computes the minimum of the condition number over all elements.

Parameters: mesh (fenics.Mesh) – The mesh, whose quality shall be computed.
Returns: The minimal condition number quality measure.
Return type: float

classmethod min_maximum_angle(mesh)¶

Computes the minimal quality measure based on the largest angle.

This measures the relative distance of a triangle’s angles or a tetrahedron’s dihedral angles to the corresponding optimal angle. The optimal angle is defined as the angle an equilateral (and thus equiangular) element has. This is defined as

\[1 - \max\left( \frac{\alpha - \alpha^*}{\pi - \alpha^*} , 0 \right),\]

where \(\alpha\) is the corresponding (dihedral) angle of the element and \(\alpha^*\) is the corresponding (dihedral) angle of the reference element. To compute the quality measure, the minimum of this expression over all elements and all of their (dihedral) angles is computed.

Parameters: mesh (fenics.Mesh) – The mesh, whose quality shall be computed.
Returns: The minimum value of the maximum angle quality measure.
Return type: float

static min_radius_ratios(mesh)¶

Computes the minimal radius ratio of the mesh.

This measures the ratio of the element’s inradius to it’s circumradius, normalized by the geometric dimension. This is computed via

\[d \frac{r}{R},\]

where \(d\) is the spatial dimension, \(r\) is the inradius, and \(R\) is the circumradius. To compute the (global) quality measure, the minimum of this expression over all elements is returned.

Parameters: mesh (fenics.Mesh) – The mesh, whose quality shall be computed.
Returns: The minimal radius ratio of the mesh.
Return type: float

classmethod min_skewness(mesh)¶

Computes the minimal skewness of the mesh.

This measure the relative distance of a triangle’s angles or a tetrahedron’s dihedral angles to the corresponding optimal angle. The optimal angle is defined as the angle an equilateral, and thus equiangular, element has. The skewness lies in \([0,1]\), where 1 corresponds to the case of an optimal (equilateral) element, and 0 corresponds to a degenerate element. The skewness corresponding to some (dihedral) angle \(\alpha\) is defined as

\[1 - \max \left( \frac{\alpha - \alpha^*}{\pi - \alpha*} , \frac{\alpha^* - \alpha}{\alpha^* - 0} \right),\]

where \(\alpha^*\) is the corresponding angle of the reference element. To compute the quality measure, the minimum of this expression over all elements and all of their (dihedral) angles is computed.

Parameters: mesh (fenics.Mesh) – The mesh whose quality shall be computed.
Returns: The skewness of the mesh.
Return type: float

Functions¶

The functions which are directly available in cashocs are taken from the sub-modules geometry, nonlinear_solvers, and utils. These are functions that are likely to be used often, so that they are directly callable via cashocs.function for any of the functions shown below. Note, that they are repeated in the API reference for their respective sub-modules.

import_mesh¶

cashocs.import_mesh(input_arg)¶

Imports a mesh file for use with CASHOCS / FEniCS.

This function imports a mesh file that was generated by GMSH and converted to .xdmf with the command line function cashocs-convert. If there are Physical quantities specified in the GMSH file, these are imported to the subdomains and boundaries output of this function and can also be directly accessed via the measures, e.g., with dx(1), ds(1), etc.

Parameters

input_arg (str or configparser.ConfigParser) – This is either a string, in which case it corresponds to the location of the mesh file in .xdmf file format, or a config file that has this path stored in its settings, under the section Mesh, as parameter mesh_file.

Returns

mesh (Mesh) – The imported (computational) mesh.
subdomains (fenics.MeshFunction) – A fenics.MeshFunction object containing the subdomains, i.e., the Physical regions marked in the GMSH file.
boundaries (fenics.MeshFunction) – A MeshFunction object containing the boundaries, i.e., the Physical regions marked in the GMSH file. Can, e.g., be used to set up boundary conditions.
dx (fenics.Measure) – The volume measure of the mesh corresponding to the subdomains (i.e. GMSH Physical region indices).
ds (fenics.Measure) – The surface measure of the mesh corresponding to the boundaries (i.e. GMSH Physical region indices).
dS (fenics.Measure) – The interior facet measure of the mesh corresponding to boundaries (i.e. GMSH Physical region indices).

Return type

Tuple[cashocs.geometry.Mesh, fenics.MeshFunction, fenics.MeshFunction, fenics.Measure, fenics.Measure, fenics.Measure]

Notes

In case the boundaries in the Gmsh .msh file are not only marked with numbers (as pyhsical groups), but also with names (i.e. strings), these strings can be used with the integration measures dx and ds returned by this method. E.g., if one specified the following in a 2D Gmsh .geo file

Physical Surface("domain", 1) = {i,j,k};

where i,j,k are representative for some integers, then this can be used in the measure dx (as we are 2D) as follows. The command

dx(1)

is completely equivalent to

dx("domain")

and both can be used interchangeably.

regular_mesh¶

cashocs.regular_mesh(n=10, L_x=1.0, L_y=1.0, L_z=None, diagonal='right')¶

Creates a mesh corresponding to a rectangle or cube.

This function creates a uniform mesh of either a rectangle or a cube, starting at the origin and having length specified in L_x, L_y, and L_z. The resulting mesh uses n elements along the shortest direction and accordingly many along the longer ones. The resulting domain is

\[\begin{split}\begin{alignedat}{2} &[0, L_x] \times [0, L_y] \quad &&\text{ in } 2D, \\ &[0, L_x] \times [0, L_y] \times [0, L_z] \quad &&\text{ in } 3D. \end{alignedat}\end{split}\]

The boundary markers are ordered as follows:

1 corresponds to \(x=0\).

2 corresponds to \(x=L_x\).

3 corresponds to \(y=0\).

4 corresponds to \(y=L_y\).

5 corresponds to \(z=0\) (only in 3D).

6 corresponds to \(z=L_z\) (only in 3D).

Parameters

n (int) – Number of elements in the shortest coordinate direction.
L_x (float) – Length in x-direction.
L_y (float) – Length in y-direction.
L_z (float or None, optional) – Length in z-direction, if this is None, then the geometry will be two-dimensional (default is None).
diagonal (str, optional) – This defines the type of diagonal used to create the box mesh in 2D. This can be one of "right", "left", "left/right", "right/left" or "crossed".

Returns

mesh (fenics.Mesh) – The computational mesh.
subdomains (fenics.MeshFunction) – A fenics.MeshFunction object containing the subdomains.
boundaries (fenics.MeshFunction) – A MeshFunction object containing the boundaries.
dx (fenics.Measure) – The volume measure of the mesh corresponding to subdomains.
ds (fenics.Measure) – The surface measure of the mesh corresponding to boundaries.
dS (fenics.Measure) – The interior facet measure of the mesh corresponding to boundaries.

Return type

Tuple[fenics.Mesh, fenics.MeshFunction, fenics.MeshFunction, fenics.Measure, fenics.Measure, fenics.Measure]

regular_box_mesh¶

cashocs.regular_box_mesh(n=10, S_x=0.0, S_y=0.0, S_z=None, E_x=1.0, E_y=1.0, E_z=None, diagonal='right')¶

Creates a mesh corresponding to a rectangle or cube.

This function creates a uniform mesh of either a rectangle or a cube, with specified start (S_) and end points (E_). The resulting mesh uses n elements along the shortest direction and accordingly many along the longer ones. The resulting domain is

\[\begin{split}\begin{alignedat}{2} &[S_x, E_x] \times [S_y, E_y] \quad &&\text{ in } 2D, \\ &[S_x, E_x] \times [S_y, E_y] \times [S_z, E_z] \quad &&\text{ in } 3D. \end{alignedat}\end{split}\]

The boundary markers are ordered as follows:

1 corresponds to \(x=S_x\).

2 corresponds to \(x=E_x\).

3 corresponds to \(y=S_y\).

4 corresponds to \(y=E_y\).

5 corresponds to \(z=S_z\) (only in 3D).

6 corresponds to \(z=E_z\) (only in 3D).

Parameters

n (int) – Number of elements in the shortest coordinate direction.
S_x (float) – Start of the x-interval.
S_y (float) – Start of the y-interval.
S_z (float or None, optional) – Start of the z-interval, mesh is 2D if this is None (default is None).
E_x (float) – End of the x-interval.
E_y (float) – End of the y-interval.
E_z (float or None, optional) – End of the z-interval, mesh is 2D if this is None (default is None).
diagonal (str, optional) – This defines the type of diagonal used to create the box mesh in 2D. This can be one of "right", "left", "left/right", "right/left" or "crossed".

Returns

mesh (fenice.Mesh) – The computational mesh.
subdomains (fenics.MeshFunction) – A MeshFunction object containing the subdomains.
boundaries (fenics.MeshFunction) – A MeshFunction object containing the boundaries.
dx (fenics.Measure) – The volume measure of the mesh corresponding to subdomains.
ds (fenics.Measure) – The surface measure of the mesh corresponding to boundaries.
dS (fenics.Measure) – The interior facet measure of the mesh corresponding to boundaries.

Return type

Tuple[fenics.Mesh, fenics.MeshFunction, fenics.MeshFunction, fenics.Measure, fenics.Measure, fenics.Measure]

load_config¶

cashocs.load_config(path)¶

Loads a config object from a config file.

Loads the config from a .ini file via the configparser package.

Parameters: path (str) – The path to the .ini file storing the configuration.
Returns: The output config file, which includes the path to the .ini file.
Return type: configparser.ConfigParser

create_dirichlet_bcs¶

cashocs.create_dirichlet_bcs(function_space, value, boundaries, idcs, **kwargs)¶

Create several Dirichlet boundary conditions at once.

Wraps multiple Dirichlet boundary conditions into a list, in case they have the same value but are to be defined for multiple boundaries with different markers. Particularly useful for defining homogeneous boundary conditions.

Parameters

function_space (fenics.FunctionSpace) – The function space onto which the BCs should be imposed on.
value (fenics.Constant or fenics.Expression or fenics.Function or float or tuple(float)) – The value of the boundary condition. Has to be compatible with the function_space, so that it could also be used as fenics.DirichletBC(function_space, value, ...).
boundaries (fenics.MeshFunction) – The fenics.MeshFunction object representing the boundaries.
idcs (list[int or str] or int or str) – A list of indices / boundary markers that determine the boundaries onto which the Dirichlet boundary conditions should be applied to. Can also be a single entry for a single boundary. If your mesh file is named, then you can also use the names of the boundaries to define the boundary conditions.
**kwargs – Keyword arguments for fenics.DirichletBC

Returns

A list of DirichletBC objects that represent the boundary conditions.

Return type

list[fenics.DirichletBC]

Examples

Generate homogeneous Dirichlet boundary conditions for all 4 sides of the unit square

from fenics import *
import cashocs

mesh, _, _, _, _, _ = cashocs.regular_mesh(25)
V = FunctionSpace(mesh, 'CG', 1)
bcs = cashocs.create_dirichlet_bcs(V, Constant(0), boundaries, [1,2,3,4])

newton_solve¶

cashocs.newton_solve(F, u, bcs, dF=None, shift=None, rtol=1e-10, atol=1e-10, max_iter=50, convergence_type='combined', norm_type='l2', damped=True, inexact=True, verbose=True, ksp=None, ksp_options=None, A_tensor=None, b_tensor=None)¶

A damped Newton method for solving nonlinear equations.

The damped Newton method is based on the natural monotonicity test from Deuflhard, Newton methods for nonlinear problems. It also allows fine tuning via a direct interface, and absolute, relative, and combined stopping criteria. Can also be used to specify the solver for the inner (linear) subproblems via petsc ksps.

The method terminates after max_iter iterations, or if a termination criterion is satisfied. These criteria are given by

a relative one in case convergence_type = 'rel', i.e.,

\[\lvert\lvert F_{k} \rvert\rvert \leq \texttt{rtol} \lvert\lvert F_0 \rvert\rvert.\]

an absolute one in case convergence_type = 'abs', i.e.,

\[\lvert\lvert F_{k} \rvert\rvert \leq \texttt{atol}.\]

a combination of both in case convergence_type = 'combined', i.e.,

\[\lvert\lvert F_{k} \rvert\rvert \leq \texttt{atol} + \texttt{rtol} \lvert\lvert F_0 \rvert\rvert.\]

The norm chosen for the termination criterion is specified via norm_type.

Parameters

F (ufl.form.Form) – The variational form of the nonlinear problem to be solved by Newton’s method.
u (fenics.Function) – The sought solution / initial guess. It is not assumed that the initial guess satisfies the Dirichlet boundary conditions, they are applied automatically. The method overwrites / updates this Function.
bcs (list[fenics.DirichletBC]) – A list of DirichletBCs for the nonlinear variational problem.
dF (ufl.Form or None, optional) – The Jacobian of F, used for the Newton method. Default is None, and in this case the Jacobian is computed automatically with AD.
shift (ufl.Form or None, optional) – This is used in case we want to solve a nonlinear operator equation with a nonlinear part F and a part shift, which does not depend on the variable u. Solves the equation \(F(u) = shift\). In case shift is None (the default), the equation \(F(u) = 0\) is solved.
rtol (float, optional) – Relative tolerance of the solver if convergence_type is either 'combined' or 'rel' (default is rtol = 1e-10).
atol (float, optional) – Absolute tolerance of the solver if convergence_type is either 'combined' or 'abs' (default is atol = 1e-10).
max_iter (int, optional) – Maximum number of iterations carried out by the method (default is max_iter = 50).
convergence_type ({'combined', 'rel', 'abs'}) – Determines the type of stopping criterion that is used.
norm_type ({'l2', 'linf'}) – Determines which norm is used in the stopping criterion.
damped (bool, optional) – If True, then a damping strategy is used. If False, the classical Newton-Raphson iteration (without damping) is used (default is True).
inexact (bool, optional) – If True, then an inexact Newtons method is used, in case an iterative solver is used for the inner solution of the linear systems. Default is True.
verbose (bool, optional) – If True, prints status of the iteration to the console (default is True).
ksp (petsc4py.PETSc.KSP or None, optional) – The PETSc ksp object used to solve the inner (linear) problem if this is None it uses the direct solver MUMPS (default is None).
ksp_options (list[list[str]] or None, optional) – The list of options for the linear solver.
A_tensor (Optional[fenics.PETScMatrix]) –
b_tensor (Optional[fenics.PETScVector]) –

Returns

The solution of the nonlinear variational problem, if converged. This overrides the input function u.

Return type

fenics.Function

Examples

Consider the problem

\[\begin{split}\begin{alignedat}{2} - \Delta u + u^3 &= 1 \quad &&\text{ in } \Omega=(0,1)^2 \\ u &= 0 \quad &&\text{ on } \Gamma. \end{alignedat}\end{split}\]

This is solved with the code

from fenics import *
import cashocs

mesh, _, boundaries, dx, _, _ = cashocs.regular_mesh(25)
V = FunctionSpace(mesh, 'CG', 1)

u = Function(V)
v = TestFunction(V)
F = inner(grad(u), grad(v))*dx + pow(u,3)*v*dx - Constant(1)*v*dx
bcs = cashocs.create_dirichlet_bcs(V, Constant(0.0), boundaries, [1,2,3,4])
cashocs.damped_newton_solve(F, u, bcs)

set_log_level¶

cashocs.set_log_level(level)¶

Determines the log level of cashocs.

Can be used to show, e.g., info and warning messages or to hide them. There are a total of five different levels for the logs: DEBUG, INFO, WARNING, ERROR, and CRITICAL. The usage of this method is explained in the examples section.

Parameters: level (int) – Should be one of cashocs.LogLevel.DEBUG, cashocs.LogLevel.INFO, cashocs.LogLevel.WARNING, cashocs.LogLevel.ERROR, cashocs.LogLevel.CRITICAL
Returns
Return type: None

Notes

The log level setting is global, so if you use this interactively, you have to restart / reload your interactive console to return to the default settings.

Examples

To set the log level of cashocs, use this method as follows:

import cashocs

cashocs.set_log_level(cashocs.LogLevel.WARNING)

which only shows messages with a level of WARNING or higher. To use a different level, replace WARNING by DEBUG, INFO, ERROR, or CRITICAL.

Deprecated Capabilities¶

Here, we list deprecated functions and classes of cashocs, which can still be used, but are replaced by newer instances with more capabilities. These deprecated objects are still maintained for compatibility reasons.

damped_newton_solve¶

cashocs.damped_newton_solve(F, u, bcs, dF=None, rtol=1e-10, atol=1e-10, max_iter=50, convergence_type='combined', norm_type='l2', damped=True, verbose=True, ksp=None, ksp_options=None)¶

Damped Newton solve interface, only here for compatibility reasons.

Parameters

F (ufl.form.Form) – The variational form of the nonlinear problem to be solved by Newton’s method.
u (fenics.Function) – The sought solution / initial guess. It is not assumed that the initial guess satisfies the Dirichlet boundary conditions, they are applied automatically. The method overwrites / updates this Function.
bcs (list[fenics.DirichletBC]) – A list of DirichletBCs for the nonlinear variational problem.
dF (ufl.form.Form, optional) – The Jacobian of F, used for the Newton method. Default is None, and in this case the Jacobian is computed automatically with AD.
rtol (float, optional) – Relative tolerance of the solver if convergence_type is either 'combined' or 'rel' (default is rtol = 1e-10).
atol (float, optional) – Absolute tolerance of the solver if convergence_type is either 'combined' or 'abs' (default is atol = 1e-10).
max_iter (int, optional) – Maximum number of iterations carried out by the method (default is max_iter = 50).
convergence_type ({'combined', 'rel', 'abs'}) – Determines the type of stopping criterion that is used.
norm_type ({'l2', 'linf'}) – Determines which norm is used in the stopping criterion.
damped (bool, optional) – If True, then a damping strategy is used. If False, the classical Newton-Raphson iteration (without damping) is used (default is True).
verbose (bool, optional) – If True, prints status of the iteration to the console (default is True).
ksp (petsc4py.PETSc.KSP, optional) – The PETSc ksp object used to solve the inner (linear) problem if this is None it uses the direct solver MUMPS (default is None).
ksp_options (list[list[str]]) – The list of options for the linear solver.

Returns

The solution of the nonlinear variational problem, if converged. This overrides the input function u.

Return type

fenics.Function

Deprecated since version 1.5.0: This is replaced by cashocs.newton_solve and will be removed in the future.

create_config¶

cashocs.create_config(path)¶

Loads a config object from a config file.

Loads the config from a .ini file via the configparser package.

Parameters: path (str) – The path to the .ini file storing the configuration.
Returns: The output config file, which includes the path to the .ini file.
Return type: configparser.ConfigParser

Deprecated since version 1.1.0: This is replaced by load_config and will be removed in the future.

create_bcs_list¶

cashocs.create_bcs_list(function_space, value, boundaries, idcs, **kwargs)¶

Create several Dirichlet boundary conditions at once.

Wraps multiple Dirichlet boundary conditions into a list, in case they have the same value but are to be defined for multiple boundaries with different markers. Particularly useful for defining homogeneous boundary conditions.

Parameters

function_space (fenics.FunctionSpace) – The function space onto which the BCs should be imposed on.
value (fenics.Constant or fenics.Expression or fenics.Function or float or tuple(float)) – The value of the boundary condition. Has to be compatible with the function_space, so that it could also be used as fenics.DirichletBC(function_space, value, ...).
boundaries (fenics.MeshFunction) – The fenics.MeshFunction object representing the boundaries.
idcs (list[int] or int) – A list of indices / boundary markers that determine the boundaries onto which the Dirichlet boundary conditions should be applied to. Can also be a single integer for a single boundary.
**kwargs – Keyword arguments for fenics.DirichletBC

Returns

A list of DirichletBC objects that represent the boundary conditions.

Return type

list[fenics.DirichletBC]

Examples

Generate homogeneous Dirichlet boundary conditions for all 4 sides of the unit square

from fenics import *
import cashocs

mesh, _, _, _, _, _ = cashocs.regular_mesh(25)
V = FunctionSpace(mesh, 'CG', 1)
bcs = cashocs.create_dirichlet_bcs(V, Constant(0), boundaries, [1,2,3,4])

Deprecated since version 1.5.0: This is replaced by cashocs.create_dirichlet_bcs and will be removed in the future.

Sub-Modules¶

CASHOCS’ sub-modules include several additional classes and methods that could be potentially useful for the user. For the corresponding API documentation, we include the previously detailed objects, too, as to give a complete documentation of the sub-module.