BilinearOperator
A bilinear operator allows to add matrices to the system matrix that usually refer to linearisations of the PDE operators or stabilisations. If the bilinear operator lives on face entities, also jumps of operators can be involved, if they are naturally continuous for the finite element space in operation (also jumps for broken spaces) and only involve degrees of freedom on the face, e.g. normal jumps for Hdiv spaces or jumps for H1-conforming spaces or tangential jumps of Hcurl spaces. For all other discontinuous operator evaluations (that needs to evaluat more than the degrees of freedom on the face) there is the possibility to use BilinearOperatorDG. It is also possible to assign a matrix assembled by the user as a BilinearOperator.
Constructors
ExtendableFEM.BilinearOperator
— Typefunction BilinearOperator(
A::AbstractMatrix,
u_test,
u_ansatz = u_test;
kwargs...)
Generates a bilinear form from a user-provided matrix, which can be a sparse matrix or a FEMatrix with multiple blocks. The arguments utest and uansatz specify where to put the (blocks of the) matrix in the system.
ExtendableFEM.BilinearOperator
— Typefunction BilinearOperator(
[kernel!::Function],
oa_test::Array{<:Tuple{Union{Unknown,Int}, DataType},1},
oa_ansatz::Array{<:Tuple{Union{Unknown,Int}, DataType},1} = oa_test;
kwargs...)
Generates a bilinear form that evaluates the vector product of the operator evaluation(s) of the test function(s) with the operator evaluation(s) of the ansatz function(s). If a function is provided in the first argument, the ansatz function evaluations can be customized by the kernel function and its result vector is then used in a dot product with the test function evaluations. In this case the header of the kernel functions needs to be conform to the interface
kernel!(result, eval_ansatz, qpinfo)
where qpinfo allows to access information at the current quadrature point.
Operator evaluations are tuples that pair an unknown identifier or integer with a Function operator.
Example: BilinearOperator([grad(1)], [grad(1)]; kwargs...) generates a weak Laplace operator.
Keyword arguments:
bonus_quadorder
: additional quadrature order added to quadorder. Default: 0entities
: assemble operator on these grid entities (default = ONCELLS). Default: ONCELLSentry_tolerance
: threshold to add entry to sparse matrix. Default: 0factor
: factor that should be multiplied during assembly. Default: 1lump
: diagonal lumping (=0 no lumping, =1 only keep diagonal entry, =2 accumulate full row to diagonal). Default: 0name
: name for operator used in printouts. Default: ''BilinearOperator''parallel
: assemble operator in parallel using colors/partitions information (assembles into full matrix directly). Default: falseparallel_groups
: assemble operator in parallel using CellAssemblyGroups (assembles seperated matrices that are added together sequantially). Default: falseparams
: array of parameters that should be made available in qpinfo argument of kernel function. Default: nothingquadorder
: quadrature order. Default: ''auto''regions
: subset of regions where operator should be assembly only. Default: Any[]store
: store matrix separately (and copy from there when reassembly is triggered). Default: falsetime_dependent
: operator is time-dependent ?. Default: falsetransposed_copy
: assemble a transposed copy of that operator into the transposed matrix block(s), 0 = no, 1 = symmetric, -1 = skew-symmetric. Default: 0use_sparsity_pattern
: read sparsity pattern of jacobian of kernel to find out which components couple. Default: ''auto''verbosity
: verbosity level. Default: 0
ExtendableFEM.BilinearOperator
— Methodfunction BilinearOperator(
kernel::Function,
oa_test::Array{<:Tuple{Union{Unknown,Int}, DataType},1},
oa_ansatz::Array{<:Tuple{Union{Unknown,Int}, DataType},1},
oa_args::Array{<:Tuple{Union{Unknown,Int}, DataType},1};
kwargs...)
Generates a nonlinear bilinear form that evaluates a kernel function that depends on the operator evaluation(s) of the ansatz function(s) and the operator evaluations of the current solution. The result of the kernel function is used in a vector product with the operator evaluation(s) of the test function(s). Hence, this can be used as a linearization of a nonlinear operator. The header of the kernel functions needs to be conform to the interface
kernel!(result, eval_ansatz, eval_args, qpinfo)
where qpinfo allows to access information at the current quadrature point.
Operator evaluations are tuples that pair an unknown identifier or integer with a Function operator.
Example: BilinearOperator([grad(1)], [grad(1)]; kwargs...) generates a weak Laplace operator.
Keyword arguments:
bonus_quadorder
: additional quadrature order added to quadorder. Default: 0entities
: assemble operator on these grid entities (default = ONCELLS). Default: ONCELLSentry_tolerance
: threshold to add entry to sparse matrix. Default: 0factor
: factor that should be multiplied during assembly. Default: 1lump
: diagonal lumping (=0 no lumping, =1 only keep diagonal entry, =2 accumulate full row to diagonal). Default: 0name
: name for operator used in printouts. Default: ''BilinearOperator''parallel
: assemble operator in parallel using colors/partitions information (assembles into full matrix directly). Default: falseparallel_groups
: assemble operator in parallel using CellAssemblyGroups (assembles seperated matrices that are added together sequantially). Default: falseparams
: array of parameters that should be made available in qpinfo argument of kernel function. Default: nothingquadorder
: quadrature order. Default: ''auto''regions
: subset of regions where operator should be assembly only. Default: Any[]store
: store matrix separately (and copy from there when reassembly is triggered). Default: falsetime_dependent
: operator is time-dependent ?. Default: falsetransposed_copy
: assemble a transposed copy of that operator into the transposed matrix block(s), 0 = no, 1 = symmetric, -1 = skew-symmetric. Default: 0use_sparsity_pattern
: read sparsity pattern of jacobian of kernel to find out which components couple. Default: ''auto''verbosity
: verbosity level. Default: 0
BilinearOperatorDG
BilinearOperatorDG is intended for bilinear forms that involves jumps of discontinuous quantities on faces whose assembly requires evaluation of all degrees of freedom on the neighbouring cells, e.g. gradient jumps for H1 conforming functions. In this case the assembly loop triggers integration along the boundary of the cells.
ExtendableFEM.BilinearOperatorDG
— Typefunction BilinearOperatorDG(
[kernel!::Function],
oa_test::Array{<:Tuple{Union{Unknown,Int}, DataType},1},
oa_ansatz::Array{<:Tuple{Union{Unknown,Int}, DataType},1} = oa_test;
kwargs...)
Generates a bilinear form that evaluates the vector product of the (discontinuous) operator evaluation(s) of the test function(s) with the (discontinuous) operator evaluation(s) of the ansatz function(s). If a function is provided in the first argument, the ansatz function evaluations can be customized by the kernel function and its result vector is then used in a dot product with the test function evaluations. In this case the header of the kernel functions needs to be conform to the interface
kernel!(result, eval_ansatz, qpinfo)
where qpinfo allows to access information at the current quadrature point.
Operator evaluations are tuples that pair an unknown identifier or integer with a Function operator.
Example: BilinearOperatorDG([jump(grad(1))], [jump(grad(1))]; kwargs...) generates an interior penalty stabilisation.
Keyword arguments:
bonus_quadorder
: additional quadrature order added to quadorder. Default: 0callback!
: function with interface (A, b, sol) that is called in each assembly step. Default: nothingentities
: assemble operator on these grid entities (default = ONFACES). Default: ONFACESentry_tolerance
: threshold to add entry to sparse matrix. Default: 0factor
: factor that should be multiplied during assembly. Default: 1lump
: lump the operator (= only assemble the diagonal). Default: falsename
: name for operator used in printouts. Default: ''BilinearOperatorDG''parallel_groups
: assemble operator in parallel using CellAssemblyGroups. Default: falseparams
: array of parameters that should be made available in qpinfo argument of kernel function. Default: nothingquadorder
: quadrature order. Default: ''auto''regions
: subset of regions where operator should be assembly only. Default: Any[]store
: store matrix separately (and copy from there when reassembly is triggered). Default: falsetime_dependent
: operator is time-dependent ?. Default: falsetransposed_copy
: assemble a transposed copy of that operator into the transposed matrix block(s), 0 = no, 1 = symmetric, -1 = skew-symmetric. Default: 0use_sparsity_pattern
: read sparsity pattern of jacobian of kernel to find out which components couple. Default: ''auto''verbosity
: verbosity level. Default: 0
ExtendableFEM.BilinearOperatorDG
— Methodfunction BilinearOperatorDG(
kernel::Function,
oa_test::Array{<:Tuple{Union{Unknown,Int}, DataType},1},
oa_ansatz::Array{<:Tuple{Union{Unknown,Int}, DataType},1},
oa_args::Array{<:Tuple{Union{Unknown,Int}, DataType},1};
kwargs...)
Generates a nonlinear bilinear form that evaluates a kernel function that depends on the (discontinuou) operator evaluation(s) of the ansatz function(s) and the (discontinuous) operator evaluations of the current solution. The result of the kernel function is used in a vector product with the operator evaluation(s) of the test function(s). Hence, this can be used as a linearization of a nonlinear operator. The header of the kernel functions needs to be conform to the interface
kernel!(result, eval_ansatz, eval_args, qpinfo)
where qpinfo allows to access information at the current quadrature point.
Operator evaluations are tuples that pair an unknown identifier or integer with a Function operator.
Keyword arguments:
bonus_quadorder
: additional quadrature order added to quadorder. Default: 0entities
: assemble operator on these grid entities (default = ONCELLS). Default: ONCELLSentry_tolerance
: threshold to add entry to sparse matrix. Default: 0factor
: factor that should be multiplied during assembly. Default: 1lump
: diagonal lumping (=0 no lumping, =1 only keep diagonal entry, =2 accumulate full row to diagonal). Default: 0name
: name for operator used in printouts. Default: ''BilinearOperator''parallel
: assemble operator in parallel using colors/partitions information (assembles into full matrix directly). Default: falseparallel_groups
: assemble operator in parallel using CellAssemblyGroups (assembles seperated matrices that are added together sequantially). Default: falseparams
: array of parameters that should be made available in qpinfo argument of kernel function. Default: nothingquadorder
: quadrature order. Default: ''auto''regions
: subset of regions where operator should be assembly only. Default: Any[]store
: store matrix separately (and copy from there when reassembly is triggered). Default: falsetime_dependent
: operator is time-dependent ?. Default: falsetransposed_copy
: assemble a transposed copy of that operator into the transposed matrix block(s), 0 = no, 1 = symmetric, -1 = skew-symmetric. Default: 0use_sparsity_pattern
: read sparsity pattern of jacobian of kernel to find out which components couple. Default: ''auto''verbosity
: verbosity level. Default: 0
Examples
Below two examples illustrate some use cases.
Example - Stokes operator
For the linear operator of a Stokes problem a kernel could look like
μ = 0.1 # viscosity parameter
function kernel!(result, input, qpinfo)
∇u, p = view(input,1:4), view(input, 5)
result[1] = μ*∇u[1] - p[1]
result[2] = μ*∇u[2]
result[3] = μ*∇u[3]
result[4] = μ*∇u[4] - p[1]
result[5] = -(∇u[1] + ∇u[4])
return nothing
end
and the coressponding BilinearOperator constructor call reads
u = Unknown("u"; name = "velocity")
p = Unknown("p"; name = "pressure")
BilinearOperator(kernel!, [grad(u), id(p)]; use_sparsity_pattern = true)
The additional argument causes that the zero pressure-pressure block of the matrix is not (even tried to be) assembled, since input[5]
does not couple with result[5]
.
Example - interior penalty stabilization
A popular convection stabilization is based on the jumps of the gradient, which can be realised with the kernel
function stab_kernel!(result, input, qpinfo)
result .= input .* qpinfo.volume^2
end
and the BilinearOperatorDG constructor call
u = Unknown("u")
assign_operator!(PD, BilinearOperatorDG(stab_kernel!, [jump(grad(u))]; entities = ON_IFACES, factor = 0.01))