GitHub

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 — Type
function 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.

source
ExtendableFEM.BilinearOperator — Type
function 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: 0

  • entities: assemble operator on these grid entities (default = ONCELLS). Default: ONCELLS

  • entry_tolerance: threshold to add entry to sparse matrix. Default: 0

  • factor: factor that should be multiplied during assembly. Default: 1

  • lump: diagonal lumping (=0 no lumping, =1 only keep diagonal entry, =2 accumulate full row to diagonal). Default: 0

  • name: name for operator used in printouts. Default: ''BilinearOperator''

  • parallel: assemble operator in parallel using colors/partitions information (assembles into full matrix directly). Default: false

  • parallel_groups: assemble operator in parallel using CellAssemblyGroups (assembles separated matrices that are added together sequantially). Default: false

  • params: array of parameters that should be made available in qpinfo argument of kernel function. Default: nothing

  • quadorder: 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: false

  • time_dependent: operator is time-dependent ?. Default: false

  • transposed_copy: assemble a transposed copy of that operator into the transposed matrix block(s), 0 = no, 1 = symmetric, -1 = skew-symmetric. Default: 0

  • use_sparsity_pattern: read sparsity pattern of jacobian of kernel to find out which components couple. Default: ''auto''

  • verbosity: verbosity level. Default: 0

source
ExtendableFEM.BilinearOperator — Method
function 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: 0

  • entities: assemble operator on these grid entities (default = ONCELLS). Default: ONCELLS

  • entry_tolerance: threshold to add entry to sparse matrix. Default: 0

  • factor: factor that should be multiplied during assembly. Default: 1

  • lump: diagonal lumping (=0 no lumping, =1 only keep diagonal entry, =2 accumulate full row to diagonal). Default: 0

  • name: name for operator used in printouts. Default: ''BilinearOperator''

  • parallel: assemble operator in parallel using colors/partitions information (assembles into full matrix directly). Default: false

  • parallel_groups: assemble operator in parallel using CellAssemblyGroups (assembles separated matrices that are added together sequantially). Default: false

  • params: array of parameters that should be made available in qpinfo argument of kernel function. Default: nothing

  • quadorder: 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: false

  • time_dependent: operator is time-dependent ?. Default: false

  • transposed_copy: assemble a transposed copy of that operator into the transposed matrix block(s), 0 = no, 1 = symmetric, -1 = skew-symmetric. Default: 0

  • use_sparsity_pattern: read sparsity pattern of jacobian of kernel to find out which components couple. Default: ''auto''

  • verbosity: verbosity level. Default: 0

source

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 — Type
function 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: 0

  • callback!: function with interface (A, b, sol) that is called in each assembly step. Default: nothing

  • entities: assemble operator on these grid entities (default = ONFACES). Default: ONFACES

  • entry_tolerance: threshold to add entry to sparse matrix. Default: 0

  • factor: factor that should be multiplied during assembly. Default: 1

  • lump: lump the operator (= only assemble the diagonal). Default: false

  • name: name for operator used in printouts. Default: ''BilinearOperatorDG''

  • parallel: assemble operator in parallel using colors/partitions information (assembles into full matrix directly). Default: false

  • parallel_groups: assemble operator in parallel using CellAssemblyGroups (assembles separated matrices that are added together sequantially). Default: false

  • params: array of parameters that should be made available in qpinfo argument of kernel function. Default: nothing

  • quadorder: 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: false

  • time_dependent: operator is time-dependent ?. Default: false

  • transposed_copy: assemble a transposed copy of that operator into the transposed matrix block(s), 0 = no, 1 = symmetric, -1 = skew-symmetric. Default: 0

  • use_sparsity_pattern: read sparsity pattern of jacobian of kernel to find out which components couple. Default: ''auto''

  • verbosity: verbosity level. Default: 0

source
ExtendableFEM.BilinearOperatorDG — Method
function 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: 0

  • entities: assemble operator on these grid entities (default = ONCELLS). Default: ONCELLS

  • entry_tolerance: threshold to add entry to sparse matrix. Default: 0

  • factor: factor that should be multiplied during assembly. Default: 1

  • lump: diagonal lumping (=0 no lumping, =1 only keep diagonal entry, =2 accumulate full row to diagonal). Default: 0

  • name: name for operator used in printouts. Default: ''BilinearOperator''

  • parallel: assemble operator in parallel using colors/partitions information (assembles into full matrix directly). Default: false

  • parallel_groups: assemble operator in parallel using CellAssemblyGroups (assembles separated matrices that are added together sequantially). Default: false

  • params: array of parameters that should be made available in qpinfo argument of kernel function. Default: nothing

  • quadorder: 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: false

  • time_dependent: operator is time-dependent ?. Default: false

  • transposed_copy: assemble a transposed copy of that operator into the transposed matrix block(s), 0 = no, 1 = symmetric, -1 = skew-symmetric. Default: 0

  • use_sparsity_pattern: read sparsity pattern of jacobian of kernel to find out which components couple. Default: ''auto''

  • verbosity: verbosity level. Default: 0

source

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))