Postprocessing

Plotting

Plotting can be performed using the package GridVisualize.jl. This package extends the API with a couple of methods:

GridVisualize.scalarplotFunction
scalarplot(
    sys::VoronoiFVM.AbstractSystem,
    sol::AbstractMatrix;
    species,
    scale,
    kwargs...
) -> Any

Plot one species from solution

source
scalarplot(
    sys::VoronoiFVM.AbstractSystem,
    sol::VoronoiFVM.AbstractTransientSolution;
    species,
    scale,
    tscale,
    kwargs...
) -> Any

Plot one species from transient solution

source
GridVisualize.scalarplot!Function
scalarplot!(
    vis,
    sys::VoronoiFVM.AbstractSystem,
    sol::AbstractMatrix;
    species,
    scale,
    kwargs...
) -> Any

Plot one species from solution

source
scalarplot!(
    vis,
    sys::VoronoiFVM.AbstractSystem,
    sol::VoronoiFVM.AbstractTransientSolution;
    species,
    scale,
    tscale,
    tlabel,
    kwargs...
) -> Any

Plot one species from transient solution

source
VoronoiFVM.plothistoryFunction
plothistory(tsol; 
            plots=[:timestepsizes,
                   :timestepupdates,
                   :newtonsteps,
                   :newtonupdates], 
            size=(700,600), 
            logmin=1.0e-20,
            Plotter=GridVisualize.default_plotter(),
            kwargs...)

Plot solution history stored in tsol. The plots argument allows to choose which plots are shown.

source

Norms & volumes

LinearAlgebra.normFunction
norm(system, u)
norm(system, u, p)

Calculate Euklidean norm of the degree of freedom vector.

source
VoronoiFVM.lpnormFunction
lpnorm(sys, u, p)
lpnorm(sys, u, p, species_weights)

Calculate weighted discrete $L^p$ norm of a solution vector.

source
VoronoiFVM.l2normFunction
l2norm(sys, u)
l2norm(sys, u, species_weights)

Calculate weigthed discrete $L^2(\Omega)$ norm of a solution vector.

source
VoronoiFVM.w1pseminormFunction
w1pseminorm(sys, u, p)
w1pseminorm(sys, u, p, species_weights)

Calculate weighted discrete $W^{1,p}(\Omega)$ seminorm of a solution vector.

source
VoronoiFVM.h1seminormFunction
h1seminorm(sys, u)
h1seminorm(sys, u, species_weights)

Calculate weighted discrete $H^1(\Omega)$ seminorm of a solution vector.

source
VoronoiFVM.w1pnormFunction
w1pnorm(sys, u, p)
w1pnorm(sys, u, p, species_weights)

Calculate weighted discrete $W^{1,p}(\Omega)$ norm of a solution vector.

source
VoronoiFVM.h1normFunction
h1norm(sys, u)
h1norm(sys, u, species_weights)

Calculate weighted discrete $H^1(\Omega)$ norm of a solution vector.

source
VoronoiFVM.lpw1pseminormFunction
lpw1pseminorm(sys, u, p)
lpw1pseminorm(sys, u, p, species_weights)

Calculate weighted discrete $L^p([0,T];W^{1,p}(\Omega))$ seminorm of a transient solution.

source
VoronoiFVM.l2h1seminormFunction
l2h1seminorm(sys, u)
l2h1seminorm(sys, u, species_weights)

Calculate weighted discrete $L^2([0,T];H^1(\Omega))$ seminorm of a transient solution.

source
VoronoiFVM.lpw1pnormFunction
lpw1pnorm(sys, u, p)
lpw1pnorm(sys, u, p, species_weights)

Calculate weighted discrete $L^p([0,T];W^{1,p}(\Omega))$ norm of a transient solution.

source
VoronoiFVM.l2h1normFunction
l2h1norm(sys, u)
l2h1norm(sys, u, species_weights)

Calculate weighted discrete $L^2([0,T];H^1(\Omega))$ norm of a transient solution.

source

Solution integrals

VoronoiFVM.integrateMethod
integrate(system,F,U; boundary=false)

Integrate node function (same signature as reaction or storage) F of solution vector region-wise over domain or boundary. The result is an nspec x nregion vector.

source
VoronoiFVM.edgeintegrateFunction
edgeintegrate(system,F,U; boundary=false)

Integrate edge function (same signature as flux function) F of solution vector region-wise over domain or boundary. The result is an nspec x nregion vector.

source

Nodal flux reconstruction

VoronoiFVM.nodefluxFunction
nodeflux(system, U)

Reconstruction of edge flux as vector function on the nodes of the triangulation. The result can be seen as a piecewiesw linear vector function in the FEM space spanned by the discretization nodes exhibiting the flux density.

The reconstruction is based on the "magic formula" R. Eymard, T. Gallouet, R. Herbin, IMA Journal of Numerical Analysis (2006) 26, 326−353, Lemma 2.4 .

The return value is a dim x nspec x nnodes vector. The flux of species i can e.g. plotted via GridVisualize.vectorplot.

Example:

    ispec=3
    vis=GridVisualizer(Plotter=Plotter)
    scalarplot!(vis,grid,solution[ispec,:],clear=true,colormap=:summer)
    vectorplot!(vis,grid,nf[:,ispec,:],clear=false)
    reveal(vis)

CAVEAT: there is a possible unsolved problem with the values at domain corners in the code. Please see any potential boundary artifacts as a manifestation of this issue and report them.

source

Boundary flux calculation

VoronoiFVM.TestFunctionFactoryType
mutable struct TestFunctionFactory

Data structure containing DenseSystem used to calculate test functions for boundary flux calculations.

  • system::VoronoiFVM.AbstractSystem: Original system
  • tfsystem::VoronoiFVM.System{Tv, Tc, Ti, Tm, Matrix{Ti}, Matrix{Tv}} where {Tv, Tc, Ti, Tm}: Test function system
  • control::SolverControl: Solver control
source
VoronoiFVM.TestFunctionFactoryMethod
TestFunctionFactory(
    system::VoronoiFVM.AbstractSystem;
    control
) -> TestFunctionFactory

Constructor for TestFunctionFactory from System

source
VoronoiFVM.integrateMethod
integrate(system, tf, U; kwargs...)

Calculate test function integral for steady state solution.

source
VoronoiFVM.integrateMethod
integrate(system, tf, U, Uold, tstep; params)

Calculate test function integral for transient solution.

source
VoronoiFVM.testfunctionMethod
testfunction(factory::TestFunctionFactory, bc0, bc1) -> Any

Create testfunction which has Dirichlet zero boundary conditions for boundary regions in bc0 and Dirichlet one boundary conditions for boundary regions in bc1.

source

Impedance calculatiom

Impedance calculation can be seen as a postprocessing step after the solution of the unexcited stationary system.

VoronoiFVM.ImpedanceSystemType
mutable struct ImpedanceSystem{Tv} <: VoronoiFVM.AbstractImpedanceSystem{Tv}

Concrete type for impedance system.

  • sysnzval::AbstractArray{Complex{Tv}, 1} where Tv: Nonzero pattern of time domain system matrix
  • storderiv::AbstractMatrix: Derivative of storage term
  • matrix::AbstractArray{Complex{Tv}, 2} where Tv: Complex matrix of impedance system
  • F::AbstractArray{Complex{Tv}, 2} where Tv: Right hand side of impedance system
  • U0::AbstractMatrix: Stationary state
source
VoronoiFVM.ImpedanceSystemMethod
ImpedanceSystem(system, U0, excited_spec, excited_bc)

Construct impedance system from time domain system sys and steady state solution U0 under the assumption of a periodic perturbation of species excited_spec at boundary excited_bc.

source
VoronoiFVM.impedanceMethod
impedance(impedance_system,ω, U0 ,
          excited_spec, excited_bc, excited_bcval,
           dmeas_stdy,
           dmeas_tran 
           )
    

Calculate impedance.

  • ω: frequency
  • U0: steady state slution
  • dmeas_stdy: Derivative of steady state part of measurement functional
  • dmeas_tran Derivative of transient part of the measurement functional
source
VoronoiFVM.measurement_derivativeMethod
measurement_derivative(system, measurement_functional, U0)

Calculate the derivative of the scalar measurement functional at steady state U0

Usually, this functional is a test function integral. Initially, we assume that its value depends on all unknowns of the system.

source