proteus.mprans.VOF module

class proteus.mprans.VOF.SubgridError(coefficients, nd)[source]

Bases: proteus.SubgridError.SGE_base

initializeElementQuadrature(mesh, t, cq)[source]
updateSubgridErrorHistory(initializationPhase=False)[source]
calculateSubgridError(q)[source]
class proteus.mprans.VOF.ShockCapturing(coefficients, nd, shockCapturingFactor=0.25, lag=True, nStepsToDelay=None)[source]

Bases: proteus.ShockCapturing.ShockCapturing_base

initializeElementQuadrature(mesh, t, cq)[source]
updateShockCapturingHistory()[source]
class proteus.mprans.VOF.NumericalFlux(vt, getPointwiseBoundaryConditions, getAdvectiveFluxBoundaryConditions, getDiffusiveFluxBoundaryConditions)[source]

Bases: proteus.NumericalFlux.Advection_DiagonalUpwind_Diffusion_IIPG_exterior

class proteus.mprans.VOF.RKEV(transport, timeOrder=1, runCFL=0.1, integrateInterpolationPoints=False)[source]

Bases: proteus.TimeIntegration.SSP

choose_dt()[source]
initialize_dt(t0, tOut, q)[source]

Modify self.dt

setCoefficients()[source]

beta are all 1’s here mwf not used right now

updateStage()[source]

Need to switch to use coefficients

initializeTimeHistory(resetFromDOF=True)[source]

Push necessary information into time history arrays

updateTimeHistory(resetFromDOF=False)[source]

assumes successful step has been taken

generateSubsteps(tList)[source]

create list of substeps over time values given in tList. These correspond to stages

resetOrder(order)[source]

initialize data structures for stage updges

setFromOptions(nOptions)[source]

allow classes to set various numerical parameters

TimeIntegration = <module 'proteus.TimeIntegration' from '/home/cekees/proteus/proteus/TimeIntegration.py'>[source]
class proteus.mprans.VOF.Coefficients(LS_model=None, V_model=0, RD_model=None, ME_model=1, EikonalSolverFlag=0, checkMass=True, epsFact=0.0, useMetrics=0.0, sc_uref=1.0, sc_beta=1.0, setParamsFunc=None, movingDomain=False, forceStrongConditions=0, STABILIZATION_TYPE=0, ENTROPY_TYPE=2, LUMPED_MASS_MATRIX=False, FCT=True, cE=1.0, uL=0.0, uR=1.0, cK=1.0)[source]

Bases: proteus.TransportCoefficients.TC_base

initializeMesh(mesh)[source]
attachModels(modelList)[source]
initializeElementQuadrature(t, cq)[source]
initializeElementBoundaryQuadrature(t, cebq, cebq_global)[source]
initializeGlobalExteriorElementBoundaryQuadrature(t, cebqe)[source]
preStep(t, firstStep=False)[source]
postStep(t, firstStep=False)[source]
updateToMovingDomain(t, c)[source]
evaluate(t, c)[source]
class EikonalSolver(levelSolverType, F, relativeTolerance=0.0, absoluteTolerance=1e-08, maxSolverIts=100, frontTolerance=0.0001, frontInitType='magnitudeOnly', bandTolerance=-0.01, eikonalVariable=0, localReconstruction=None, printInfo=False)[source]

Simple wrapper for special purpose Eikonal equation solvers on a single level. Current types allowed:

FMMEikonalSolver
FSWEikonalSolver
convertFromC0P1Rep(dofin, dofout)[source]

allow output FEM space to be something besides C0P1 and convert from C0P1 to expected representation here

convertToC0P1Rep(dofin, dofout, interpType='min')[source]

allow input FEM space to be something besides C0P1 and then convert to expected representation here

info()[source]
solve(u, r, b=None)[source]
class Coefficients.FMMEikonalSolver(mesh, dofMap, nSpace, localSolverType='QianEtalV2', frontInitType='magnitudeOnly', debugLevel=3)[source]

Encapsulate naive implementation of Fast Marching Methods on unstructured grids for

\[\|\grad T\| = 1/F\]

\(T = 0\) on \(\Gamma\)

1d local solver is standard upwind approximation

2d local solver variations: acute triangulations version 1 or version 2 from Qian Zhang etal 07 obtuse triangulation not implemented

3d local solver varitions: not fully checked

For now, the input should be non-negative!

cfmmfsw = <module 'proteus.cfmmfsw' from '/home/cekees/proteus/proteus/cfmmfsw.so'>[source]
solve(phi0, T, nodalSpeeds=None, zeroTol=0.0001, trialTol=0.1, verbose=0)[source]

Test first order fast marching method algorithm for eikonal equation

|grad T | = 1, phi(

ec x) = 0, x in Gamma

assuming phi_0 describes initial location of interface Gamma and has reasonable values (absolute values) for T close to Gamma. Here T can be interpreted as the travel time from Gamma.

Right now assumes global node numbers <–> global dofs but this can be fixed easily Input

phi0: dof array from P1 C0 FiniteElementFunction holding initial condition

T : dof array from P1 C0 FiniteElementFunction for solution

Output T(

ec x_n) : travel time from initial front to node ( ec x_n)

Internal data structures

Status
: status of nodal point (dictionary)
-1 –> Far
0 –> Trial 1 –> Known

Trial : nodal points adjacent to front tuples (index,val) stored in heap

TODO
have return flag
class Coefficients.FSWEikonalSolver(mesh, dofMap, nSpace, iterAtol=1e-08, iterRtol=0.0, maxIts=100, localSolverType='QianEtalV2', frontInitType='magnitudeOnly', refPoints=None, orderApprox=1, LARGE=1.234e+28, debugLevel=3)[source]
Encapsulate naive implementation of Fast Marching Methods on unstructured grids
for
\[\|\grad T\| = 1/F\]

\(T = 0\) on \(\Gamma\)

1d local solver is standard upwind approximation

2d local solver variations: acute triangulations version 1 or version 2 from Qian Zhang etal 07 obtuse triangulation not implemented

3d local solver variations: not fully checked

For now, the input should be non-negative!

cfmmfsw = <module 'proteus.cfmmfsw' from '/home/cekees/proteus/proteus/cfmmfsw.so'>[source]
solve(phi0, T, nodalSpeeds=None, zeroTol=0.0001, trialTol=0.1, verbose=0)[source]

Test first order fast sweeping method algorithm for eikonal equation

|grad T | = 1, phi(

ec x) = 0, x in Gamma

assuming phi_0 describes initial location of interface Gamma and has reasonable values (absolute values) for T close to Gamma. Here T can be interpreted as the travel time from Gamma.

Right now assumes global node numbers <–> global dofs but this can be fixed easily Input

phi0: dof array holding P1 C0 FiniteElementFunction holding initial condition

T : dof array holding P1 C0 FiniteElementFunction for solution

Output T(

ec x_n) : travel time from initial front to node ( ec x_n)

Internal data structures

Status
: status of nodal point (dictionary)
0 –> Not Known (Trial) 1 –> Known

Order : ordering of points in domain using l_p metric from fixed reference points

Coefficients.VOFCoefficientsEvaluate()[source]

evaluate the coefficients of the conservative level set advection equation

Coefficients.VolumeAveragedVOFCoefficientsEvaluate()[source]

evaluate the coefficients of the conservative level set advection equation with variable porosity

Coefficients.copyExteriorElementBoundaryValuesFromElementBoundaryValues()[source]

copy quantity in element boundary quadrature array to one that lives only on exterior element boundaries

class proteus.mprans.VOF.LevelModel(uDict, phiDict, testSpaceDict, matType, dofBoundaryConditionsDict, dofBoundaryConditionsSetterDict, coefficients, elementQuadrature, elementBoundaryQuadrature, fluxBoundaryConditionsDict=None, advectiveFluxBoundaryConditionsSetterDict=None, diffusiveFluxBoundaryConditionsSetterDictDict=None, stressTraceBoundaryConditionsSetterDict=None, stabilization=None, shockCapturing=None, conservativeFluxDict=None, numericalFluxType=None, TimeIntegrationClass=None, massLumping=False, reactionLumping=False, options=None, name='defaultName', reuse_trial_and_test_quadrature=True, sd=True, movingDomain=False, bdyNullSpace=False)[source]

Bases: proteus.Transport.OneLevelTransport

nCalls = 0[source]
FCTStep()[source]
calculateCoefficients()[source]
updateVelocityFieldAsFunction()[source]
calculateElementResidual()[source]
getResidual(u, r)[source]
getJacobian(jacobian)[source]
calculateElementQuadrature()[source]

Calculate the physical location and weights of the quadrature rules and the shape information at the quadrature points.

This function should be called only when the mesh changes.

calculateElementBoundaryQuadrature()[source]
calculateExteriorElementBoundaryQuadrature()[source]

Calculate the physical location and weights of the quadrature rules and the shape information at the quadrature points on global element boundaries.

This function should be called only when the mesh changes.

estimate_mt()[source]
calculateSolutionAtQuadrature()[source]
calculateAuxiliaryQuantitiesAfterStep()[source]
updateAfterMeshMotion()[source]