proteus.TimeIntegration module

A class hierarchy of time discretizations

class proteus.TimeIntegration.TI_base(transport, integrateInterpolationPoints=False)

Bases: object

The base class for time discretizations of tightly coupled systems of transport equations.

The purpose of time integration objects is to construct approximations:

m_t, dm_t – the time derivatives and their derivatives w.r.t. u u_*, du_* – the solution variables and their derivatives w.r.t. u f_*, df_*,… –the coefficients and their derivatives w.r.t. u

The objects provide a set_dt function for setting the time step as well as functions for SUGGESTING the time step. The StepController objects have the final responsibility for choosing the time step (by calling set_dt).

The objects can also enforce a tolerance on an a posteriori error estimate using the lastStepErrorOK member function.

If any of the coefficient approximations are not implicit, a flag can be set to pass that information to the transport object.

This base class implements no time integration method and results in the steady state scalar transport equation:

deld ( f - a grad phi) + r = 0

NoIntegration is an alias for this class.

Set flags that indicate that all terms are implicit.

calculateU(u)

Generate u_*

calculateElementCoefficients(q)

Calculate m_t, and dm_t and recalculate any of the other coefficients on element interiors

calculateElementBoundaryCoefficients(ebq)

Calculate m_t, and dm_t and recalculate any of the other coefficients on element boundaries

calculateExteriorElementBoundaryCoefficients(ebqe)

Calculate m_t, and dm_t and recalculate any of the other coefficients on global element boundaries

calculateStrongElementSpatialResidual(q)

Recalculate the part of the element residual due to spatial terms

calculateElementSpatialResidual(elementResidual)

Recalculate the part of the element residual due to spatial terms

calculateElementSpatialJacobian(elementJacobian)

Recalculate the part of the element Jacobian due to the spatial terms

calculateElementSpatialBoundaryJacobian(elementBoundaryJacobian, elementBoundaryJacobian_eb)

Recalculate the part of the element Jacobian due to the spatial terms along element boundaries

calculateExteriorElementSpatialBoundaryJacobian(exteriorElementBoundaryJacobian)

Recalculate the part of the element Jacobian due to the spatial terms along element boundaries

calculateGeneralizedInterpolationCoefficients(cip)

Calculate m_t, and dm_t and recalculate any of the other coefficients on element interiors

initialize_dt(t0, tOut, q)

Modify self.dt

choose_dt()

Modify self.dt

set_dt(DTSET)

generateSubsteps(tList)

create list of substeps over time values given in tList.

initializeSpaceHistory(resetFromDOF=False)

updateTimeHistory(resetFromDOF=False)

Push necessary information into time history arrays

resetTimeHistory(resetFromDOF=False)

initializeTimeHistory(resetFromDOF=True)

Push necessary information into time history arrays

updateStage()

increment stage counters and push necessary information into stage arrays

setInitialStageValues()

set the stage values assuming this is the first step after a problem reset

setInitialGuess()

set an initial guess for the next step

setLastSolveFailed(lastSolveFailed)

tell integrator last attempted time step failed or not

lastStepErrorOk()

was the last time step acceptable or not

calculateCoefs()

calculate any coefficients that depend on new time step

setFromOptions(nOptions)

allow classes to set various numerical parameters

postAdaptUpdate(oldTime)

ensure array histories are properly set after adapt. This should only work for backwards euler and vbdf at the moment.

proteus.TimeIntegration.NoIntegration

alias of proteus.TimeIntegration.TI_base

class proteus.TimeIntegration.BackwardEuler(transport, integrateInterpolationPoints=False)

Bases: proteus.TimeIntegration.TI_base

Set flags that indicate that all terms are implicit.

calculateElementCoefficients(q)

Calculate m_t, and dm_t and recalculate any of the other coefficients on element interiors

calculateGeneralizedInterpolationCoefficients(cip)

Calculate m_t, and dm_t and recalculate any of the other coefficients on element interiors

initializeTimeHistory(resetFromDOF=True)

Push necessary information into time history arrays

updateTimeHistory(resetFromDOF=False)

Push necessary information into time history arrays

calculateCoefs()

calculate any coefficients that depend on new time step

class proteus.TimeIntegration.BackwardEuler_cfl(transport, runCFL=0.9, integrateInterpolationPoints=False)

Bases: proteus.TimeIntegration.BackwardEuler

Take a fraction of the max (over all components and elements)

Set flags that indicate that all terms are implicit.

choose_dt()

Modify self.dt

initialize_dt(t0, tOut, q)

Modify self.dt

updateTimeHistory(resetFromDOF=False)

Push necessary information into time history arrays

setFromOptions(nOptions)

allow classes to set various numerical parameters

postAdaptUpdate(oldTime)

This looks exactly like updateTimeHistory with the exception of setting self.tLast=self.t

class proteus.TimeIntegration.SSP(transport, runCFL=0.9, integrateInterpolationPoints=False)

Bases: proteus.TimeIntegration.BackwardEuler_cfl

Set flags that indicate that all terms are implicit.

class proteus.TimeIntegration.FLCBDF(transport)

Bases: proteus.TimeIntegration.TI_base

Set flags that indicate that all terms are implicit.

initialize_dt(t0, tOut, q)

Modify self.dt

initializeTimeHistory(resetFromDOF=False)

Push necessary information into time history arrays

calculateCoefs()

calculate any coefficients that depend on new time step

calculateElementCoefficients(q)

Calculate m_t, and dm_t and recalculate any of the other coefficients on element interiors

choose_dt()

Modify self.dt

set_dt(DTSET)

updateTimeHistory(resetFromDOF=False)

Push necessary information into time history arrays

updateStage()

increment stage counters and push necessary information into stage arrays

setInitialStageValues()

set the all stage values assuming this is the first step after a problem reset

lastStepErrorOk()

was the last time step acceptable or not

class proteus.TimeIntegration.PsiTCtte(transport)

Bases: proteus.TimeIntegration.BackwardEuler_cfl

Set flags that indicate that all terms are implicit.

psiTCtteDT()

initializeTimeHistory(resetFromDOF=True)

Push necessary information into time history arrays

calculateElementCoefficients(q)

Calculate m_t, and dm_t and recalculate any of the other coefficients on element interiors

choose_dt()

Modify self.dt

updateTimeHistory(resetFromDOF=False)

Push necessary information into time history arrays

setFromOptions(nOptions)

allow classes to set various numerical parameters

class proteus.TimeIntegration.PsiTCtte_new(transport)

Bases: proteus.TimeIntegration.BackwardEuler

Set flags that indicate that all terms are implicit.

psiTCtteDT()

calculateElementCoefficients(q)

Calculate m_t, and dm_t and recalculate any of the other coefficients on element interiors

choose_dt()

Modify self.dt

resetTimeHistory(resetFromDOF=False)

updateTimeHistory(resetFromDOF=False)

Push necessary information into time history arrays

setFromOptions(nOptions)

allow classes to set various numerical parameters

class proteus.TimeIntegration.ForwardEuler(transport, runCFL=0.45)

Bases: proteus.TimeIntegration.TI_base

Set flags that indicate that all terms are implicit.

calculateElementSpatialJacobian(elementJacobian)

Recalculate the part of the element Jacobian due to the spatial terms

calculateElementBoundaryCoefficients(ebq)

Calculate m_t, and dm_t and recalculate any of the other coefficients on element boundaries do I need these to be synchronized with q or do the elementBoundary quad values get lagged in the spatial residual

calculateElementSpatialBoundaryJacobian(elementBoundaryJacobian, elementBoundaryJacobian_eb)

Recalculate the part of the element Jacobian due to the spatial terms along element boundaries

calculateExteriorElementSpatialBoundaryJacobian(elementBoundaryJacobian)

Recalculate the part of the element Jacobian due to the spatial terms along element boundaries

calculateElementCoefficients(q)

Calculate m_t, and dm_t and recalculate any of the other coefficients on element interiors

calculateStrongElementSpatialResidual(q)

Recalculate the part of the element residual due to spatial terms

calculateElementSpatialResidual(elementResidual)

Recalculate the part of the element residual due to spatial terms

choose_dt()

Modify self.dt

initialize_dt(t0, tOut, q)

Modify self.dt

initializeSpaceHistory(resetFromDOF=False)

updateTimeHistory(resetFromDOF=False)

Push necessary information into time history arrays

setFromOptions(nOptions)

allow classes to set various numerical parameters

class proteus.TimeIntegration.ForwardEuler_A(transport, runCFL=0.45, limiterType=None)

Bases: proteus.TimeIntegration.TI_base

Set flags that indicate that all terms are implicit.

calculateElementCoefficients(q)

Calculate m_t, and dm_t and recalculate any of the other coefficients on element interiors

calculateElementBoundaryCoefficients(ebq)

Calculate m_t, and dm_t and recalculate any of the other coefficients on element boundaries

calculateExteriorElementBoundaryCoefficients(ebqe)

Calculate m_t, and dm_t and recalculate any of the other coefficients on global element boundaries

initialize_dt(t0, tOut, q)

Modify self.dt

choose_dt()

Modify self.dt mwf needs to be checked

updateStage()

if using limiting

updateTimeHistory(resetFromDOF=False)

Push necessary information into time history arrays

setFromOptions(nOptions)

allow classes to set various numerical parameters

initializeSpaceHistory(resetFromDOF=True)

initializeTimeHistory(resetFromDOF=True)

Push necessary information into time history arrays

class proteus.TimeIntegration.ForwardEuler_H(transport, runCFL=0.45)

Bases: proteus.TimeIntegration.TI_base

Set flags that indicate that all terms are implicit.

calculateU(u)

Generate u_*

calculateElementCoefficients(q)

Calculate m_t, and dm_t and recalculate any of the other coefficients on element interiors

calculateElementBoundaryCoefficients(ebq)

Calculate m_t, and dm_t and recalculate any of the other coefficients on element boundaries

calculateElementSpatialResidual(elementResidual)

Recalculate the part of the element residual due to spatial terms

choose_dt()

Modify self.dt mwf needs to be checked

updateTimeHistory(resetFromDOF=False)

Push necessary information into time history arrays

initializeSpaceHistory(resetFromDOF=False)

updateStage()

increment stage counters and push necessary information into stage arrays

setInitialStageValues()

set the all stage values assuming this is the first step after a problem reset

setFromOptions(nOptions)

allow classes to set various numerical parameters

class proteus.TimeIntegration.OuterTheta(transport, runCFL=0.9, thetaAdvection=0.5, thetaDiffusion=0.5, thetaReaction=0.5, thetaHamiltonian=0.5)

Bases: proteus.TimeIntegration.TI_base

Set flags that indicate that all terms are implicit.

calculateU(u)

Generate u_*

calculateElementCoefficients(q)

Calculate m_t, and dm_t and recalculate any of the other coefficients on element interiors

calculateElementBoundaryCoefficients(ebq)

Calculate m_t, and dm_t and recalculate any of the other coefficients on element boundaries

choose_dt()

Modify self.dt mwf needs to be checked

updateTimeHistory(resetFromDOF=False)

Push necessary information into time history arrays

initializeSpaceHistory(resetFromDOF=False)

setFromOptions(nOptions)

allow classes to set various numerical parameters

proteus.TimeIntegration.IMEX_AD_Euler

alias of proteus.TimeIntegration.ForwardEuler_A

proteus.TimeIntegration.IMEX_HJ_Euler

alias of proteus.TimeIntegration.ForwardEuler_H

class proteus.TimeIntegration.VBDF(transport, timeOrder=2, integrateInterpolationPoints=False)

Bases: proteus.TimeIntegration.TI_base

Variable coefficient bdf’s just first and second order

Set flags that indicate that all terms are implicit.

calculateElementCoefficients(q)

Calculate m_t, and dm_t and recalculate any of the other coefficients on element interiors

initialize_dt(t0, tOut, q)

This routine is bad right now, because we don’t want to compute an internal dt anymore just the error estimates

Controller now has to take care of self.t and self.dt

choose_dt()

This routine is bad right now, because we don’t want to compute an internal dt anymore just the error estimates

Controller now has to take care of self.t and self.dt

setInitialGuess()

set an initial guess for the next step need to figure out how to synchronize this

predictor necessary for error estimates

lastStepErrorOk()

was the last time step acceptable or not This routine now just computes the error estimate …

initializeTimeHistory(resetFromDOF=True)

Push necessary information into time history arrays

updateTimeHistory(resetFromDOF=False)

Push necessary information into time history arrays

postAdaptUpdate(oldTime)

This looks exactly like updateTimeHistory with the exeption of setting self.tLast=self.t

computeErrorEstimate()

calculate :math:` ec{e}^{n+1}`

depending on order of approximation To be consistent, this must be called after step taken but before update time history

Initially, use

\[(\]

ec y^{n+1}- ec y^{n+1,p})/2

for first order and

ec y^{n+1}-(1+r) ec y^{n}+r ec y^{n-1})

for second order.

calculatePredictor()

calculateCoefs()

calculate any coefficients that depend on new time step

chooseOrder()

setFromOptions(nOptions)

allow classes to set various numerical parameters

class proteus.TimeIntegration.ExplicitRK_base(transport, timeOrder=1, nStages=1, runCFL=0.1, usingDTinMass=True)

Bases: proteus.TimeIntegration.TI_base

base class for explicit RK schemes which can be written in the Shu Osher 88 (Cockburn Shu 89) style generic formula

u^0 = u^n u^l = sum_{k=0}^{l-1}lpha_{l,k}u^k + eta_{l,k} Delta t L(u^{k},t^n + d_{k} Delta t^n) u^{n+1} = u^{m}

We will try to include these schemes in a Rothe framework rather than MOL

so in calculate element coefficients, at stage number l, u_t will be replaced with

u_t –> (u - sum_{k=0}^{l-1} lpha_{l,k} u^k)/Delta t

or

u_t –> (u - sum_{k=0}^{l-1} lpha_{l,k} u^k)

depending on whether Delta t is assumed to be in the mass matrix or not

The spatial residual at stage l will be replaced with (in calculateElementSpatialResidual)

r –> sum_{k=0,l-1} eta_{l,k} r_k

or

r –> Delta t sum_{k=0,l-1} eta_{l,k} r_k

depending on whether Delta t is assumed to be in the mass matrix or not Here r_k corresponds to

L(u^{k},t^n + d_{k} Delta t^n)

Substeps in NumericalSolution framework correspond to stage values in the RK integration. These are setup in generateSubsteps. The indexing is a little confusing, since the residuals are ‘cached’ when they come into calculateElementResidual in order to be used at the next stage. This means the time levels in substep have to be for the following stage level. This is modulo the number of stages so t^{n+1} at the last stage –> t^n for the next time step

Set flags that indicate that all terms are implicit.

set_dt(DTSET)

resetOrder(order, stages)

resize and reinitialize if necessary

setCoefficients()

set alpha,beta, d for

\[u^0 = u^n u^l = \sum_{k=0}^{l-1}lpha_{l,k}u^k + eta_{l,k} \Delta t L(u^{k},t^n + d_{k} \Delta t^n) u^{n+1} = u^{m}\]

must be implemented in derived class

Note that

\[alpha[l,k] --> alpha_{l+1,k} beta[l,k] --> beta_{l+1,k} dcoefs[k] --> d_{k+1}\]

because of the whole cachingdelayed eval deal second index (k) is actual level

updateStage()

increment stage counter by 1. lstage here is last stage completed

calculateElementCoefficients(q)

set m_t as

m_t –> (m - sum_{k=0}^{l-1} lpha_{l,k} m^k)/Delta t

calculateElementCoefficientsNoDtInMass(q)

set m_t as

m_t –> (m - sum_{k=0}^{l-1} lpha_{l,k} m^k)

calculateElementSpatialResidual(elementResidual)

The spatial residual at stage l will be replaced with (in calculateElementSpatialResidual)

r –> sum_{k=0,l-1} eta_{l,k} r_k

calculateElementSpatialResidualNoDtInMass(elementResidual)

The spatial residual at stage l will be replaced with (in calculateElementSpatialResidual)

r –> Delta t sum_{k=0,l-1} eta_{l,k} r_k

calculateElementSpatialJacobian(elementJacobian)

Recalculate the part of the element Jacobian due to the spatial terms

calculateElementBoundaryCoefficients(ebq)

Calculate m_t, and dm_t and recalculate any of the other coefficients on element boundaries do I need these to be synchronized with q or do the elementBoundary quad values get lagged in the spatial residual

calculateElementSpatialBoundaryJacobian(elementBoundaryJacobian, elementBoundaryJacobian_eb)

Recalculate the part of the element Jacobian due to the spatial terms along element boundaries

calculateExteriorElementSpatialBoundaryJacobian(elementBoundaryJacobian)

Recalculate the part of the element Jacobian due to the spatial terms along element boundaries

initialize_dt(t0, tOut, q)

Modify self.dt

choose_dt()

Modify self.dt

generateSubsteps(tList)

create list of substeps over time values given in tList. These correspond to stages so need to have the right number set. For the SSP RK schemes here, t^m = t^{n+1} for all the stages except the first, I believe

The order here is a little confusing. Since the residuals are ‘cached’ when they come into calculateElementResidual in order to be used at the next stage, the time levels have to be for the next stage This is modulo the number of stages so t^{n+1} at the last stage –> t^n for the next time step

setInitialStageValues()

setup all the stage values if we assume that the initial condition has been loaded into lstage+1

initializeTimeHistory(resetFromDOF=True)

Push necessary information into time history arrays

updateTimeHistory(resetFromDOF=False)

setup for next time step, cycle U^s –> U^0

setFromOptions(nOptions)

allow classes to set various numerical parameters

class proteus.TimeIntegration.LinearSSPRKintegration(transport, order=1, runCFL=0.1, usingSSPRKNewton=False)

Bases: proteus.TimeIntegration.ExplicitRK_base

implementation of Explicit RK schemes for ODE that are SSP for linear operators (aka) spatial discretizations in

\[\od{\]

ec y}{t} = mat{L} ec y

See Gottlieb, Shu, Tadmor siam review article and notes

Their formulation for this scheme is

\[ \begin{align}\begin{aligned}u^0 = u^n u^l = u^{l-1} + \Delta t L(u^{l-1}) l = 1,...,m-1 u^m = \sum_{k=0}^{l-1}lpha_{m,k}u^k + lpha_{m,m-1} \Delta t L(u^{m-1}) u^{n+1} = u^m\\so lpha_{l,l-1} = 1, eta_{l,l-1} = 1.0 l < m, eta_{m,m-1} = lph_{m,m-1}\end{aligned}\end{align} \]

Apparently,

\[d_{l} = l \Delta t l < m, d_m = 1.0\]

which doesn’t make a lot of sense for time dependent problems

Set flags that indicate that all terms are implicit.

setCoefficients()

alpha matrix in Gottlieb, Shu, Tadmor review paper

u^{l} = sum_{k=0}^{l-1}lpha_{l,k}u^k

so lpha_{l,l-1} = 1 for l < m

beta matrix multiplies

eta_{l,k} Delta t L(u^{k},t^n + d_{k} Delta t^n)

so

eta_{l,l-1} = 1.0 for l < m eta_{m,m-1} = lpha_{m,m-1}

Note that

alpha[l,k] –> alpha_{l+1,k} beta[l,k] –> beta_{l+1,k} dcoefs[k] –> d_{k+1} because of the whole cachingdelayed eval deal

second index (k) is actual level

setFromOptions(nOptions)

allow classes to set various numerical parameters

class proteus.TimeIntegration.SSPRKPIintegration(transport, order=1, runCFL=0.1, limiterType=None, usingSSPRKNewton=False)

Bases: proteus.TimeIntegration.ExplicitRK_base

So called SSP RK integration with limiting step at end of each

phase This one holds the original TVD RK schemes up to 3rd order 1st order is Forward Euler

2nd order

\[u^1 = u^n + \Delta t L(u^n) u^2 = 0.5(u^n + u^1) + 0.5\Delta t L(u^1,t^n + \Delta t^n), same stage values etc linear SSPRK\]

3rd order

\[u^1 = u^n + \Delta t L(u^n,t^n) u^2 = 0.75 u^n + 0.25u^1 + 0.25\Delta t L(u^1,t^n + \Delta t^n), u^3 = 1/3 u^n + 2/3 u^2 + 2/3 \Delta t L(u^2,t^n + 0.5\Delta t^n)\]

generic formula is

\[u^0 = u^n u^l = \sum_{k=0}^{l-1}lpha_{l,k}u^k + eta_l \Delta t L(u^{l-1},t^n + d_{l-1} \Delta t^n) u^{n+1} = u^{m}\]

so \(beta_l --> beta_{l,k-1}\) in our basic framework

Try to put this into Rothe paradigm with linear implicit mass matrix evaluation by saying at stage l

\[m_t pprox\]

rac{u^l - sum_{k=0}^{l-1}lpha_{l,k}u^k}{Delta t}

spatial residual –> spatial residual * beta_l

and solve system for \(u^l\)

The limiting right now assumes same space for each component

Set flags that indicate that all terms are implicit.

setCoefficients()

See Cockburn and Hou and Shu

For TVD schemes beta_{l,k} = 0 k < l-1

Note that

alpha[l,k] –> alpha_{l+1,k} beta[l,k] –> beta_{l+1,k} dcoefs[k] –> d_{k+1} because of the whole cachingdelayed eval deal

updateStage()

increment stage counter by 1. lstage here is last stage completed

setFromOptions(nOptions)

allow classes to set various numerical parameters

class proteus.TimeIntegration.LinearSSPRKPIintegration(transport, order=1, runCFL=0.1, limiterType=None, usingSSPRKNewton=False)

Bases: proteus.TimeIntegration.LinearSSPRKintegration

Linear SSP RK integration with limiting step at end of each phase, this one just uses linear SSP RK coefficients so not really correct NL SSP (TVD) RK for 3rd order (or higher)

Way limiter right now assumes same space for everything

Later on, pass in limiter object I guess

updateStage()

increment stage counter by 1. lstage here is last stage completed

setFromOptions(nOptions)

allow classes to set various numerical parameters

class proteus.TimeIntegration.DGlimiterP1Lagrange1d(mesh, nSpace, u, transport=None, limiterFlag=0)

Bases: object

canonical (I hope) 1d DG limiting procedure when original local approximation space has P1 Lagrange basis assumes all solution fem spaces are identical

applyDGlimitingP1Lagrange1d()

initializeMeshInfo()

setFromOptions(nOptions)

applySlopeLimiting(uIn, uDofOut)

Apply limiting procedure directly using dofs

class proteus.TimeIntegration.DGlimiterP2Lagrange1d(mesh, nSpace, u, transport=None, limiterFlag=0)

Bases: object

canonical (I hope) 1d DG limiting procedure when original local approximation space has P2 Lagrange basis

applyDGlimitingP1Lagrange1d()

initializeMeshInfo()

setFromOptions(nOptions)

projectToLimitedSpace(solndofs, ci=0)

project fem function from femSpaceSolution with dofs held in solndofs to ulim which is DG_AffineQuadraticOnSimplexWithNodalBasis

nodal interpolant will not preserve mass necessarily

projectFromLimitedSpace(solndofs, limdofs, tag, ci=0)

project to fem function from femSpaceSolution with dofs held in solndofs from ulim which is DG_AffineQuadraticOnSimplexWithNodalBasis Don’t overwrite solution if tag == 1

start with nodal interpolant, since this will preserve mass

applySlopeLimiting(uIn, uDofOut)

Apply limiting procedure directly using dofs

class proteus.TimeIntegration.DGlimiterPkMonomial1d(mesh, nSpace, u, transport=None, limiterFlag=0)

Bases: object

canonical (I hope) 1d DG limiting procedure when original local approximation space has Pk monomial basis

applyDGlimitingP1Lagrange1d()

initializeMeshInfo()

setFromOptions(nOptions)

projectToLimitedSpace(solndofs, ci=0)

project fem function from femSpaceSolution with dofs held in solndofs to ulim which is DG_LinearQuadraticOnSimplexWithNodalBasis

try more or less generic L2 projection on reference element

projectFromLimitedSpace(solndofs, limdofs, tag, ci)

project to fem function from femSpaceSolution with dofs held in solndofs from ulim which is DG_AffineLinearOnSimplexWithNodalBasis Don’t overwrite solution if tag == 1

applySlopeLimiting(uIn, uDofOut)

Apply limiting procedure directly using dofs

class proteus.TimeIntegration.UnstructuredLimiter_base(mesh, nSpace)

Bases: object

computeElementNeighborShapeGradients()

initializeMeshInfo(verbose=0)

setFromOptions(nOptions)

class proteus.TimeIntegration.CockburnNotesLimiter2d_base(mesh, nSpace)

Bases: proteus.TimeIntegration.UnstructuredLimiter_base

computeCockburnDGlimiterArrays2d()

computeAlphaCoefs(verbose=0)

loop through each element and element face midpoint, compute barycentric coordinates for face midpoint using triangles formed by neighboring barycenters accept first nonnegative pair of coordinates

write lambda_i = 1 + (x-x_i).grad lambda_i

class proteus.TimeIntegration.DGlimiterP1Lagrange2d(mesh, nSpace, u, transport=None, nu=1.5, M=0.0)

Bases: proteus.TimeIntegration.CockburnNotesLimiter2d_base

apply standard, Cockburn RKDG limiter in 2d (maybe)

applyCockburnDGlimiterP1Lagrange2d()

applySlopeLimiting(uIn, uDofOut)

setFromOptions(nOptions)

class proteus.TimeIntegration.DGlimiterP2Lagrange2d(mesh, nSpace, u, transport=None, nu=1.5, M=0.0)

Bases: proteus.TimeIntegration.CockburnNotesLimiter2d_base

canonical (I hope) 2d DG limiting procedure when original local approximation space has P2 Lagrange basis go ahead and use a generic L2 projection though

TODO

move projection steps to c

applyCockburnDGlimiterP1Lagrange2d()

projectToLimitedSpace(solndofs, ci=0)

project fem function from femSpaceSolution with dofs held in solndofs to ulim which is DG_AffineLinearOnSimplexWithNodalBasis

try more or less generic L2 projection on reference element

projectFromLimitedSpace(solndofs, limdofs, tag, ci=0)

project to fem function from femSpaceSolution with dofs held in solndofs from ulim which is DG_AffineLinearOnSimplexWithNodalBasis Don’t overwrite solution if tag == 1

applySlopeLimiting(uIn, uDofOut)

Apply limiting procedure directly using dofs

setFromOptions(nOptions)

class proteus.TimeIntegration.DGlimiterPkMonomial2d(mesh, nSpace, u, transport=None, nu=1.5, M=0.0)

Bases: proteus.TimeIntegration.CockburnNotesLimiter2d_base

canonical (I hope) 2d DG limiting procedure when original local approximation space has Pk monomial basis

applyCockburnDGlimiterP1Lagrange2d()

projectToLimitedSpace(solndofs, ci=0)

project fem function from femSpaceSolution with dofs held in solndofs to ulim which is DG_AffineQuadraticOnSimplexWithNodalBasis

try more or less generic L2 projection on reference element

projectFromLimitedSpace(solndofs, limdofs, tag, ci=0)

project to fem function from femSpaceSolution with dofs held in solndofs from ulim which is DG_AffineQuadraticOnSimplexWithNodalBasis Don’t overwrite solution if tag == 1

start with nodal interpolant, since this will preserve mass

applySlopeLimiting(uIn, uDofOut)

Apply limiting procedure directly using dofs

setFromOptions(nOptions)

class proteus.TimeIntegration.DGlimiterDurlofskyP1Lagrange2d(mesh, nSpace, u, transport=None, killExtrema=1, allowMinWithUndershoot=0)

Bases: proteus.TimeIntegration.UnstructuredLimiter_base

applyDurlofskyDGlimiterP1Lagrange2d()

applySlopeLimiting(uIn, uDofOut)

setFromOptions(nOptions)

class proteus.TimeIntegration.DGlimiterDurlofskyP1Lagrange3d(mesh, nSpace, u, transport=None, killExtrema=1, allowMinWithUndershoot=1)

Bases: proteus.TimeIntegration.UnstructuredLimiter_base

applyDurlofskyDGlimiterP1Lagrange3d()

applySlopeLimiting(uIn, uDofOut)

setFromOptions(nOptions)

class proteus.TimeIntegration.DGlimiterP1Lagrange1d_Sw(mesh, nSpace, u, transport=None, limiterFlag=0, h_eps=0.001)

Bases: proteus.TimeIntegration.DGlimiterP1Lagrange1d

DGP1 version of limiting that applies minmod limiting for h < h_eps

and more aggressive limiting elsewhere

applyDGlimitingP1Lagrange1d_withVacuumTol()

applySlopeLimiting(uIn, uDofOut)

Apply limiting procedure directly using dofs using standard approach Then go through and eliminate negative values of water height

for cells that have average (h_bar < h_eps)

if average height is negative, then zero

if both vertices are positive leave alone (could kill slope)

if one of the vertices is negative choose slope so that this value is exactly zero

zero discharge at this vertex

May need to add additional step that limits discharge where h is much less than h_eps

setFromOptions(nOptions)

class proteus.TimeIntegration.DGlimiterP2Lagrange1d_Sw(mesh, nSpace, u, transport=None, limiterFlag=0, h_eps=0.01)

Bases: proteus.TimeIntegration.DGlimiterP2Lagrange1d

DGP1 version of limiting that applies minmod limiting for h < h_eps

and more aggressive limiting elsewhere

applyDGlimitingP1Lagrange1d_withVacuumTol()

applySlopeLimiting(uIn, uDofOut)

Apply limiting procedure directly using dofs using standard approach Then go through and eliminate negative values of water height

for cells that have average (h_bar < h_eps)

if average height is negative, then zero

if both vertices are positive leave alone (could kill slope)

if one of the vertices is negative, choose slope so that this value is exactly zero, zero discharge at this vertex

May need to add additional step that limits discharge where h is much less than h_eps

setFromOptions(nOptions)

class proteus.TimeIntegration.DGlimiterPkMonomial1d_Sw(mesh, nSpace, u, transport=None, limiterFlag=0, h_eps=0.1)

Bases: proteus.TimeIntegration.DGlimiterPkMonomial1d

DGPk version of limiting that applies minmod limiting for h < h_eps

and more aggressive limiting elsewhere

applyDGlimitingP1Lagrange1d_withVacuumTol()

applySlopeLimiting(uIn, uDofOut)

Apply limiting procedure directly using dofs

setFromOptions(nOptions)

class proteus.TimeIntegration.DGlimiterDurlofskyP1Lagrange2d_Sw(mesh, nSpace, u, transport=None, killExtrema=1, allowMinWithUndershoot=0, h_eps=0.02)

Bases: proteus.TimeIntegration.DGlimiterDurlofskyP1Lagrange2d

applyDurlofskyDGlimiterP1Lagrange2d_withVacuumTol()

applySlopeLimiting(uIn, uDofOut)

setFromOptions(nOptions)

class proteus.TimeIntegration.ForwardIntegrator(mlvtran, mlnl, dtMeth, nOptions, stepExact=False)

Bases: object

class that is responsible for basic process of integrating a problem forward in time given a VectorTranport Problem, Nonlinear Solver, and a Time Integration method

mlvTran — multilevel vector transport object for system being integrated mlNL — multilevel nonlinear solver to solve discrete system dtMet — time integration method to use nOptions — configuration options

initialize(DTSET=None, t0=0.0, T=1.0)

calculateSolution(tIn, tOut)

Move forward from time tIn to time tOut For now doesn’t worry about potential mismatch between tIn and last time value used by model

writeProgress(tn, dt, T)

just echo to screen what new and final time levels are

class proteus.TimeIntegration.SteadyStateIntegrator(mlvtran, mlnl, dtMeth, nOptions, ssatol=0.001, ssrtol=0.001, maxsteps=50, stepExact=True, ignoreFailure=False)

Bases: proteus.TimeIntegration.ForwardIntegrator

apply some form of pseudo-transient continuation or direct solve to compute steady state solution regardless of requested time step

mlvTran — multilevel vector transport object for system being integrated mlNL — multilevel nonlinear solver to solve discrete system dtMet — time integration method to use nOptions — configuration options

initialize(DTSET=None, t0=0.0, T=1.0)

calculateSolution(tIn, tOut)

integrate to steady state regardless of what tIn and tOut are

writeProgress(tn, dt, res)

just echo to screen what new and final time levels are

class proteus.TimeIntegration.SignedDistanceIntegrator(mlvtran, mlnl, dtMeth, nOptions, stepExact=False, nSteps=3)

Bases: proteus.TimeIntegration.ForwardIntegrator

apply some form of pseudo-transient continuation or direct solve to compute steady state solution regardless of requested time step

mlvTran — multilevel vector transport object for system being integrated mlNL — multilevel nonlinear solver to solve discrete system dtMet — time integration method to use nOptions — configuration options

initialize(DTSET=None, t0=0.0, T=1.0)

calculateSolution(tIn, tOut)

take a couple of steps toward steady state

writeProgress(tn, dt, res)

just echo to screen what new and final time levels are

class proteus.TimeIntegration.FLCBDF_TwophaseDarcy_fc(transport)

Bases: proteus.TimeIntegration.FLCBDF

Set flags that indicate that all terms are implicit.

calculateElementCoefficients(q)

Calculate m_t, and dm_t and recalculate any of the other coefficients on element interiors

class proteus.TimeIntegration.CentralDifference_2ndD(transport, integrateInterpolationPoints=False)

Bases: proteus.TimeIntegration.TI_base

Set flags that indicate that all terms are implicit.

calculateElementCoefficients(q)

Calculate m_t, and dm_t and recalculate any of the other coefficients on element interiors

calculateGeneralizedInterpolationCoefficients(cip)

Calculate m_t, and dm_t and recalculate any of the other coefficients on element interiors

initializeTimeHistory(resetFromDOF=True)

Push necessary information into time history arrays

updateTimeHistory(resetFromDOF=False)

Push necessary information into time history arrays

calculateCoefs()

calculate any coefficients that depend on new time step