pymor.discretizers.builtin.cg
¶
This module provides some operators for continuous finite element discretizations.
Module Contents¶
- class pymor.discretizers.builtin.cg.AdvectionOperatorP1(grid, boundary_info, advection_function=None, advection_constant=None, dirichlet_clear_columns=False, dirichlet_clear_diag=False, solver_options=None, name=None)[source]¶
Bases:
pymor.operators.numpy.NumpyMatrixBasedOperator
Linear advection
Operator
for linear finite elements.The operator is of the form
(Lu)(x) = c ∇ ⋅ [ v(x) u(x) ]
The function
v
has to be vector-valued.- Parameters:
grid – The
Grid
for which to assemble the operator.boundary_info –
BoundaryInfo
for the treatment of Dirichlet boundary conditions.advection_function – The
Function
v(x)
withshape_range = (grid.dim, )
. IfNone
, constant one is assumed.advection_constant – The constant
c
. IfNone
,c
is set to one.dirichlet_clear_columns – If
True
, set columns of the system matrix corresponding to Dirichlet boundary DOFs to zero to obtain a symmetric system matrix. Otherwise, only the rows will be set to zero.dirichlet_clear_diag – If
True
, also set diagonal entries corresponding to Dirichlet boundary DOFs to zero. Otherwise they are set to one.solver_options – The
solver_options
for the operator.name – Name of the operator.
- class pymor.discretizers.builtin.cg.AdvectionOperatorQ1(grid, boundary_info, advection_function=None, advection_constant=None, dirichlet_clear_columns=False, dirichlet_clear_diag=False, solver_options=None, name=None)[source]¶
Bases:
pymor.operators.numpy.NumpyMatrixBasedOperator
Linear advection
Operator
for bilinear finite elements.The operator is of the form
(Lu)(x) = c ∇ ⋅ [ v(x) u(x) ]
The function
v
has to be vector-valued.- Parameters:
grid – The
Grid
for which to assemble the operator.boundary_info –
BoundaryInfo
for the treatment of Dirichlet boundary conditions.advection_function – The
Function
v(x)
withshape_range = (grid.dim, )
. IfNone
, constant one is assumed.advection_constant – The constant
c
. IfNone
,c
is set to one.dirichlet_clear_columns – If
True
, set columns of the system matrix corresponding to Dirichlet boundary DOFs to zero to obtain a symmetric system matrix. Otherwise, only the rows will be set to zero.dirichlet_clear_diag – If
True
, also set diagonal entries corresponding to Dirichlet boundary DOFs to zero. Otherwise they are set to one.solver_options – The
solver_options
for the operator.name – Name of the operator.
- class pymor.discretizers.builtin.cg.BoundaryDirichletFunctional(grid, dirichlet_data, boundary_info, name=None)[source]¶
Bases:
pymor.operators.numpy.NumpyMatrixBasedOperator
Linear functional for enforcing Dirichlet boundary values.
- Parameters:
grid –
Grid
for which to assemble the functional.dirichlet_data –
Function
providing the Dirichlet boundary values.boundary_info –
BoundaryInfo
determining the Dirichlet boundaries.name – The name of the functional.
- class pymor.discretizers.builtin.cg.BoundaryL2ProductFunctional(grid, function, boundary_type=None, dirichlet_clear_dofs=False, boundary_info=None, name=None)[source]¶
Bases:
pymor.operators.numpy.NumpyMatrixBasedOperator
Linear functional representing the inner product with an L2-
Function
on the boundary.- Parameters:
grid –
Grid
for which to assemble the functional.function – The
Function
with which to take the inner product.boundary_type – The type of domain boundary (e.g. ‘neumann’) on which to assemble the functional. If
None
the functional is assembled over the whole boundary.dirichlet_clear_dofs – If
True
, set dirichlet boundary DOFs to zero.boundary_info – If
boundary_type
is specified ordirichlet_clear_dofs
isTrue
, theBoundaryInfo
determining which boundary entity belongs to which physical boundary.name – The name of the functional.
- class pymor.discretizers.builtin.cg.DiffusionOperatorP1(grid, boundary_info, diffusion_function=None, diffusion_constant=None, dirichlet_clear_columns=False, dirichlet_clear_diag=False, solver_options=None, name=None)[source]¶
Bases:
pymor.operators.numpy.NumpyMatrixBasedOperator
Diffusion
Operator
for linear finite elements.The operator is of the form
(Lu)(x) = c ∇ ⋅ [ d(x) ∇ u(x) ]
The function
d
can be scalar- or matrix-valued.- Parameters:
grid – The
Grid
for which to assemble the operator.boundary_info –
BoundaryInfo
for the treatment of Dirichlet boundary conditions.diffusion_function – The
Function
d(x)
withshape_range == ()
orshape_range = (grid.dim, grid.dim)
. IfNone
, constant one is assumed.diffusion_constant – The constant
c
. IfNone
,c
is set to one.dirichlet_clear_columns – If
True
, set columns of the system matrix corresponding to Dirichlet boundary DOFs to zero to obtain a symmetric system matrix. Otherwise, only the rows will be set to zero.dirichlet_clear_diag – If
True
, also set diagonal entries corresponding to Dirichlet boundary DOFs to zero. Otherwise they are set to one.solver_options – The
solver_options
for the operator.name – Name of the operator.
- class pymor.discretizers.builtin.cg.DiffusionOperatorQ1(grid, boundary_info, diffusion_function=None, diffusion_constant=None, dirichlet_clear_columns=False, dirichlet_clear_diag=False, solver_options=None, name=None)[source]¶
Bases:
pymor.operators.numpy.NumpyMatrixBasedOperator
Diffusion
Operator
for bilinear finite elements.The operator is of the form
(Lu)(x) = c ∇ ⋅ [ d(x) ∇ u(x) ]
The function
d
can be scalar- or matrix-valued.- Parameters:
grid – The
Grid
for which to assemble the operator.boundary_info –
BoundaryInfo
for the treatment of Dirichlet boundary conditions.diffusion_function – The
Function
d(x)
withshape_range == ()
orshape_range = (grid.dim, grid.dim)
. IfNone
, constant one is assumed.diffusion_constant – The constant
c
. IfNone
,c
is set to one.dirichlet_clear_columns – If
True
, set columns of the system matrix corresponding to Dirichlet boundary DOFs to zero to obtain a symmetric system matrix. Otherwise, only the rows will be set to zero.dirichlet_clear_diag – If
True
, also set diagonal entries corresponding to Dirichlet boundary DOFs to zero. Otherwise they are set to one.solver_options – The
solver_options
for the operator.name – Name of the operator.
- class pymor.discretizers.builtin.cg.InterpolationOperator(grid, function)[source]¶
Bases:
pymor.operators.numpy.NumpyMatrixBasedOperator
Vector-like Lagrange interpolation
Operator
for continuous finite element spaces.
- class pymor.discretizers.builtin.cg.L2ProductFunctionalP1(grid, function, dirichlet_clear_dofs=False, boundary_info=None, name=None)[source]¶
Bases:
pymor.operators.numpy.NumpyMatrixBasedOperator
Linear functional representing the inner product with an L2-
Function
.- Parameters:
grid –
Grid
for which to assemble the functional.function – The
Function
with which to take the inner product.dirichlet_clear_dofs – If
True
, set dirichlet boundary DOFs to zero.boundary_info –
BoundaryInfo
determining the Dirichlet boundaries in casedirichlet_clear_dofs
is set toTrue
.name – The name of the functional.
- class pymor.discretizers.builtin.cg.L2ProductFunctionalQ1(grid, function, dirichlet_clear_dofs=False, boundary_info=None, name=None)[source]¶
Bases:
pymor.operators.numpy.NumpyMatrixBasedOperator
Bilinear functional representing the inner product with an L2-
Function
.- Parameters:
grid –
Grid
for which to assemble the functional.function – The
Function
with which to take the inner product.dirichlet_clear_dofs – If
True
, set dirichlet boundary DOFs to zero.boundary_info –
BoundaryInfo
determining the Dirichlet boundaries in casedirichlet_clear_dofs
is set toTrue
.name – The name of the functional.
- class pymor.discretizers.builtin.cg.L2ProductP1(grid, boundary_info, dirichlet_clear_rows=True, dirichlet_clear_columns=False, dirichlet_clear_diag=False, coefficient_function=None, solver_options=None, name=None)[source]¶
Bases:
pymor.operators.numpy.NumpyMatrixBasedOperator
Operator
representing the L2-product between linear finite element functions.- Parameters:
grid – The
Grid
for which to assemble the product.boundary_info –
BoundaryInfo
for the treatment of Dirichlet boundary conditions.dirichlet_clear_rows – If
True
, set the rows of the system matrix corresponding to Dirichlet boundary DOFs to zero.dirichlet_clear_columns – If
True
, set columns of the system matrix corresponding to Dirichlet boundary DOFs to zero.dirichlet_clear_diag – If
True
, also set diagonal entries corresponding to Dirichlet boundary DOFs to zero. Otherwise, if eitherdirichlet_clear_rows
ordirichlet_clear_columns
isTrue
, the diagonal entries are set to one.coefficient_function – Coefficient
Function
for product withshape_range == ()
. IfNone
, constant one is assumed.solver_options – The
solver_options
for the operator.name – The name of the product.
- class pymor.discretizers.builtin.cg.L2ProductQ1(grid, boundary_info, dirichlet_clear_rows=True, dirichlet_clear_columns=False, dirichlet_clear_diag=False, coefficient_function=None, solver_options=None, name=None)[source]¶
Bases:
pymor.operators.numpy.NumpyMatrixBasedOperator
Operator
representing the L2-product between bilinear finite element functions.- Parameters:
grid – The
Grid
for which to assemble the product.boundary_info –
BoundaryInfo
for the treatment of Dirichlet boundary conditions.dirichlet_clear_rows – If
True
, set the rows of the system matrix corresponding to Dirichlet boundary DOFs to zero.dirichlet_clear_columns – If
True
, set columns of the system matrix corresponding to Dirichlet boundary DOFs to zero.dirichlet_clear_diag – If
True
, also set diagonal entries corresponding to Dirichlet boundary DOFs to zero. Otherwise, if eitherdirichlet_clear_rows
ordirichlet_clear_columns
isTrue
, the diagonal entries are set to one.coefficient_function – Coefficient
Function
for product withshape_range == ()
. IfNone
, constant one is assumed.solver_options – The
solver_options
for the operator.name – The name of the product.
- class pymor.discretizers.builtin.cg.RobinBoundaryOperator(grid, boundary_info, robin_data=None, solver_options=None, name=None)[source]¶
Bases:
pymor.operators.numpy.NumpyMatrixBasedOperator
Robin boundary
Operator
for linear finite elements.The operator represents the contribution of Robin boundary conditions to the stiffness matrix, where the boundary condition is supposed to be given in the form
-[ d(x) ∇u(x) ] ⋅ n(x) = c(x) (u(x) - g(x))
d
andn
are the diffusion function (seeDiffusionOperatorP1
) and the unit outer normal inx
, whilec
is the (scalar) Robin parameter function andg
is the (also scalar) Robin boundary value function.- Parameters:
grid – The
Grid
over which to assemble the operator.boundary_info –
BoundaryInfo
for the treatment of Dirichlet boundary conditions.robin_data – Tuple (c(x), g(x)) providing two
Functions
that represent the Robin parameter and boundary value function. IfNone
, the resulting operator is zero.solver_options – The
solver_options
for the operator.name – Name of the operator.
- pymor.discretizers.builtin.cg.discretize_instationary_cg(analytical_problem, diameter=None, domain_discretizer=None, grid_type=None, grid=None, boundary_info=None, num_values=None, time_stepper=None, nt=None, preassemble=True)[source]¶
Finite Element discretization of an
InstationaryProblem
.Discretizes an
InstationaryProblem
with aStationaryProblem
as the stationary part using finite elements.- Parameters:
analytical_problem – The
InstationaryProblem
to discretize.diameter – If not
None
,diameter
is passed as an argument to thedomain_discretizer
.domain_discretizer – Discretizer to be used for discretizing the analytical domain. This has to be a function
domain_discretizer(domain_description, diameter, ...)
. IfNone
,discretize_domain_default
is used.grid_type – If not
None
, this parameter is forwarded todomain_discretizer
to specify the type of the generatedGrid
.grid – Instead of using a domain discretizer, the
Grid
can also be passed directly using this parameter.boundary_info – A
BoundaryInfo
specifying the boundary types of the grid boundary entities. Must be provided ifgrid
is specified.num_values – The number of returned vectors of the solution trajectory. If
None
, each intermediate vector that is calculated is returned.time_stepper – The
time-stepper
to be used bysolve
ofInstationaryModel
.nt – If
time_stepper
is not specified, the number of time steps for implicit Euler time stepping.preassemble – If
True
, preassemble all operators in the resultingModel
.
- Returns:
m – The
Model
that has been generated.data – Dictionary with the following entries:
- grid:
The generated
Grid
.- boundary_info:
The generated
BoundaryInfo
.- unassembled_m:
In case
preassemble
isTrue
, the generatedModel
before preassembling operators.
- pymor.discretizers.builtin.cg.discretize_stationary_cg(analytical_problem, diameter=None, domain_discretizer=None, grid_type=None, grid=None, boundary_info=None, preassemble=True, mu_energy_product=None)[source]¶
Discretizes a
StationaryProblem
using finite elements.- Parameters:
analytical_problem – The
StationaryProblem
to discretize.diameter – If not
None
,diameter
is passed as an argument to thedomain_discretizer
.domain_discretizer – Discretizer to be used for discretizing the analytical domain. This has to be a function
domain_discretizer(domain_description, diameter, ...)
. IfNone
,discretize_domain_default
is used.grid_type – If not
None
, this parameter is forwarded todomain_discretizer
to specify the type of the generatedGrid
.grid – Instead of using a domain discretizer, the
Grid
can also be passed directly using this parameter.boundary_info – A
BoundaryInfo
specifying the boundary types of the grid boundary entities. Must be provided ifgrid
is specified.preassemble – If
True
, preassemble all operators in the resultingModel
.mu_energy_product – If not
None
,parameter values
for which to assemble the symmetric part of theOperator
of the resultingModel
fom
(ignoring the advection part). Thus, assuming no advection and a symmetric diffusion tensor,fom.products['energy']
is equal tofom.operator.assemble(mu)
, except for the fact that the former has cleared Dirichlet rows and columns, while the latter only has cleared Dirichlet rows).
- Returns:
m – The
Model
that has been generated.data – Dictionary with the following entries:
- grid:
The generated
Grid
.- boundary_info:
The generated
BoundaryInfo
.- unassembled_m:
In case
preassemble
isTrue
, the generatedModel
before preassembling operators.