Model Derived Type

type, public, abstract :: Model


Inherited by

type~~model~~InheritedByGraph type~model Model type~dgmodel2d_t DGModel2D_t type~dgmodel2d_t->type~model type~dgmodel1d_t DGModel1D_t type~dgmodel1d_t->type~model type~dgmodel3d_t DGModel3D_t type~dgmodel3d_t->type~model type~dgmodel3d DGModel3D type~dgmodel3d->type~dgmodel3d_t type~dgmodel2d DGModel2D type~dgmodel2d->type~dgmodel2d_t type~dgmodel2d~2 DGModel2D type~dgmodel2d~2->type~dgmodel2d_t type~dgmodel1d DGModel1D type~dgmodel1d->type~dgmodel1d_t type~dgmodel1d~2 DGModel1D type~dgmodel1d~2->type~dgmodel1d_t type~dgmodel3d~2 DGModel3D type~dgmodel3d~2->type~dgmodel3d_t type~advection_diffusion_2d_t advection_diffusion_2d_t type~advection_diffusion_2d_t->type~dgmodel2d type~linearshallowwater2d_t LinearShallowWater2D_t type~linearshallowwater2d_t->type~dgmodel2d type~nulldgmodel2d_t NullDGModel2D_t type~nulldgmodel2d_t->type~dgmodel2d type~lineareuler2d_t LinearEuler2D_t type~lineareuler2d_t->type~dgmodel2d type~nulldgmodel1d_t NullDGModel1D_t type~nulldgmodel1d_t->type~dgmodel1d type~burgers1d_t Burgers1D_t type~burgers1d_t->type~dgmodel1d type~nulldgmodel3d_t NullDGModel3D_t type~nulldgmodel3d_t->type~dgmodel3d type~advection_diffusion_3d_t advection_diffusion_3d_t type~advection_diffusion_3d_t->type~dgmodel3d type~advection_diffusion_1d_t advection_diffusion_1d_t type~advection_diffusion_1d_t->type~dgmodel1d type~linearshallowwater2d LinearShallowWater2D type~linearshallowwater2d->type~linearshallowwater2d_t type~advection_diffusion_2d~2 advection_diffusion_2d type~advection_diffusion_2d~2->type~advection_diffusion_2d_t type~lineareuler2d LinearEuler2D type~lineareuler2d->type~lineareuler2d_t type~nulldgmodel1d NullDGModel1D type~nulldgmodel1d->type~nulldgmodel1d_t type~nulldgmodel1d~2 NullDGModel1D type~nulldgmodel1d~2->type~nulldgmodel1d_t type~advection_diffusion_3d advection_diffusion_3d type~advection_diffusion_3d->type~advection_diffusion_3d_t type~advection_diffusion_3d~2 advection_diffusion_3d type~advection_diffusion_3d~2->type~advection_diffusion_3d_t type~advection_diffusion_2d advection_diffusion_2d type~advection_diffusion_2d->type~advection_diffusion_2d_t type~nulldgmodel2d NullDGModel2D type~nulldgmodel2d->type~nulldgmodel2d_t type~linearshallowwater2d~2 LinearShallowWater2D type~linearshallowwater2d~2->type~linearshallowwater2d_t type~nulldgmodel2d~2 NullDGModel2D type~nulldgmodel2d~2->type~nulldgmodel2d_t type~lineareuler2d~2 LinearEuler2D type~lineareuler2d~2->type~lineareuler2d_t type~burgers1d Burgers1D type~burgers1d->type~burgers1d_t type~burgers1d~2 Burgers1D type~burgers1d~2->type~burgers1d_t type~nulldgmodel3d NullDGModel3D type~nulldgmodel3d->type~nulldgmodel3d_t type~nulldgmodel3d~2 NullDGModel3D type~nulldgmodel3d~2->type~nulldgmodel3d_t type~advection_diffusion_1d advection_diffusion_1d type~advection_diffusion_1d->type~advection_diffusion_1d_t type~advection_diffusion_1d~2 advection_diffusion_1d type~advection_diffusion_1d~2->type~advection_diffusion_1d_t

Contents

Source Code


Components

TypeVisibilityAttributesNameInitial
real(kind=prec), public :: dt
real(kind=prec), public :: entropy
logical, public :: gradient_enabled =.false.
integer, public :: ioIterate =0
integer, public :: nvar
logical, public :: prescribed_bcs_enabled =.true.
real(kind=prec), public :: t
logical, public :: tecplot_enabled =.true.
procedure(SELF_timeIntegrator), public, pointer:: timeIntegrator=> Euler_timeIntegrator

Type-Bound Procedures

procedure, public :: AdditionalFree => AdditionalFree_Model

  • public subroutine AdditionalFree_Model(this)

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(inout) :: this

procedure, public :: AdditionalInit => AdditionalInit_Model

  • public subroutine AdditionalInit_Model(this)

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(inout) :: this

procedure, public :: AdditionalOutput => AdditionalOutput_Model

  • public subroutine AdditionalOutput_Model(this, fileid)

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(inout) :: this
    integer(kind=HID_T), intent(in) :: fileid

procedure, public :: CalculateEntropy => CalculateEntropy_Model

  • public subroutine CalculateEntropy_Model(this)

    Base method for calculating entropy of a model When this method is not overridden, the entropy is simply set to 0.0. When you develop a model built on top of this abstract class or one of its children, it is recommended that you define a convex mathematical entropy function that is used as a measure of the model stability.

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(inout) :: this

procedure(CalculateTendency), public, deferred :: CalculateTendency

  • subroutine CalculateTendency(this)Prototype

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(inout) :: this

procedure, public :: Euler_timeIntegrator

  • public subroutine Euler_timeIntegrator(this, tn)

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(inout) :: this
    real(kind=prec), intent(in) :: tn

procedure, public :: ForwardStep => ForwardStep_Model

  • public subroutine ForwardStep_Model(this, tn, dt, ioInterval)

    Forward steps the model using the associated tendency procedure and time integrator

    Read more…

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(inout) :: this
    real(kind=prec), intent(in) :: tn
    real(kind=prec), intent(in) :: dt
    real(kind=prec), intent(in) :: ioInterval

procedure, public :: GetSimulationTime

  • public subroutine GetSimulationTime(this, t)

    Returns the current simulation time stored in the model % t attribute

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(in) :: this
    real(kind=prec), intent(out) :: t

procedure, public :: IncrementIOCounter

  • public subroutine IncrementIOCounter(this)

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(inout) :: this

procedure, public :: LowStorageRK2_timeIntegrator

  • public subroutine LowStorageRK2_timeIntegrator(this, tn)

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(inout) :: this
    real(kind=prec), intent(in) :: tn

procedure, public :: LowStorageRK3_timeIntegrator

  • public subroutine LowStorageRK3_timeIntegrator(this, tn)

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(inout) :: this
    real(kind=prec), intent(in) :: tn

procedure, public :: LowStorageRK4_timeIntegrator

  • public subroutine LowStorageRK4_timeIntegrator(this, tn)

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(inout) :: this
    real(kind=prec), intent(in) :: tn

procedure, public :: PreTendency => PreTendency_Model

  • public subroutine PreTendency_Model(this)

    PreTendency is a template routine that is used to house any additional calculations that you want to execute at the beginning of the tendency calculation routine. This default PreTendency simply returns back to the caller without executing any instructions

    Read more…

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(inout) :: this

procedure, public :: PrintType => PrintType_Model

  • public subroutine PrintType_Model(this)

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(in) :: this

procedure(ReadModel), public, deferred :: ReadModel

  • subroutine ReadModel(this, filename)Prototype

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(inout) :: this
    character, intent(in) :: filename

procedure, public :: ReportEntropy => ReportEntropy_Model

  • public subroutine ReportEntropy_Model(this)

    Base method for reporting the entropy of a model to stdout. Only override this procedure if additional reporting is needed. Alternatively, if you think additional reporting would be valuable for all models, open a pull request with modifications to this base method.

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(in) :: this

procedure, public :: ReportMetrics => ReportMetrics_Model

  • public subroutine ReportMetrics_Model(this)

    Method that can be overridden by users to report their own custom metrics after file io

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(inout) :: this

procedure, public :: ReportUserMetrics => ReportUserMetrics_Model

  • public subroutine ReportUserMetrics_Model(this)

    Method that can be overridden by users to report their own custom metrics after file io

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(inout) :: this

procedure, public :: SetNumberOfVariables => SetNumberOfVariables_Model

procedure, public :: SetSimulationTime

  • public subroutine SetSimulationTime(this, t)

    Sets the model % t attribute with the provided simulation time

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(inout) :: this
    real(kind=prec), intent(in) :: t

generic, public :: SetTimeIntegrator => SetTimeIntegrator_withChar

  • public subroutine SetTimeIntegrator_withChar(this, integrator)

    Sets the time integrator method, using a character input

    Read more…

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(inout) :: this
    character, intent(in) :: integrator

procedure, private :: SetTimeIntegrator_withChar

  • public subroutine SetTimeIntegrator_withChar(this, integrator)

    Sets the time integrator method, using a character input

    Read more…

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(inout) :: this
    character, intent(in) :: integrator

procedure(UpdateGRK), public, deferred :: UpdateGRK2

  • subroutine UpdateGRK(this, m)Prototype

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(inout) :: this
    integer, intent(in) :: m

procedure(UpdateGRK), public, deferred :: UpdateGRK3

  • subroutine UpdateGRK(this, m)Prototype

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(inout) :: this
    integer, intent(in) :: m

procedure(UpdateGRK), public, deferred :: UpdateGRK4

  • subroutine UpdateGRK(this, m)Prototype

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(inout) :: this
    integer, intent(in) :: m

procedure(UpdateSolution), public, deferred :: UpdateSolution

  • subroutine UpdateSolution(this, dt)Prototype

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(inout) :: this
    real(kind=prec), intent(in), optional :: dt

procedure(WriteModel), public, deferred :: WriteModel

  • subroutine WriteModel(this, filename)Prototype

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(inout) :: this
    character, intent(in), optional :: filename

procedure(WriteTecplot), public, deferred :: WriteTecplot

  • subroutine WriteTecplot(this, filename)Prototype

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(inout) :: this
    character, intent(in), optional :: filename

procedure, public :: entropy_func => entropy_func_Model

  • public pure function entropy_func_Model(this, s) result(e)

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(in) :: this
    real(kind=prec), intent(in) :: s(1:this%nvar)

    Return Value real(kind=prec)

procedure, public :: flux1D => flux1d_Model

  • public pure function flux1d_Model(this, s, dsdx) result(flux)

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(in) :: this
    real(kind=prec), intent(in) :: s(1:this%nvar)
    real(kind=prec), intent(in) :: dsdx(1:this%nvar)

    Return Value real(kind=prec)(1:this%nvar)

procedure, public :: flux2D => flux2d_Model

  • public pure function flux2d_Model(this, s, dsdx) result(flux)

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(in) :: this
    real(kind=prec), intent(in) :: s(1:this%nvar)
    real(kind=prec), intent(in) :: dsdx(1:this%nvar,1:2)

    Return Value real(kind=prec)(1:this%nvar,1:2)

procedure, public :: flux3D => flux3d_Model

  • public pure function flux3d_Model(this, s, dsdx) result(flux)

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(in) :: this
    real(kind=prec), intent(in) :: s(1:this%nvar)
    real(kind=prec), intent(in) :: dsdx(1:this%nvar,1:3)

    Return Value real(kind=prec)(1:this%nvar,1:3)

procedure, public :: hbc1d_NoNormalFlow => hbc1d_Generic_Model

  • public pure function hbc1d_Generic_Model(this, s, nhat) result(exts)

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(in) :: this
    real(kind=prec), intent(in) :: s(1:this%nvar)
    real(kind=prec), intent(in) :: nhat

    Return Value real(kind=prec)(1:this%nvar)

procedure, public :: hbc1d_Prescribed => hbc1d_Prescribed_Model

  • public pure function hbc1d_Prescribed_Model(this, x, t) result(exts)

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(in) :: this
    real(kind=prec), intent(in) :: x
    real(kind=prec), intent(in) :: t

    Return Value real(kind=prec)(1:this%nvar)

procedure, public :: hbc1d_Radiation => hbc1d_Generic_Model

  • public pure function hbc1d_Generic_Model(this, s, nhat) result(exts)

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(in) :: this
    real(kind=prec), intent(in) :: s(1:this%nvar)
    real(kind=prec), intent(in) :: nhat

    Return Value real(kind=prec)(1:this%nvar)

procedure, public :: hbc2d_NoNormalFlow => hbc2d_Generic_Model

  • public pure function hbc2d_Generic_Model(this, s, nhat) result(exts)

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(in) :: this
    real(kind=prec), intent(in) :: s(1:this%nvar)
    real(kind=prec), intent(in) :: nhat(1:2)

    Return Value real(kind=prec)(1:this%nvar)

procedure, public :: hbc2d_Prescribed => hbc2d_Prescribed_Model

  • public pure function hbc2d_Prescribed_Model(this, x, t) result(exts)

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(in) :: this
    real(kind=prec), intent(in) :: x(1:2)
    real(kind=prec), intent(in) :: t

    Return Value real(kind=prec)(1:this%nvar)

procedure, public :: hbc2d_Radiation => hbc2d_Generic_Model

  • public pure function hbc2d_Generic_Model(this, s, nhat) result(exts)

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(in) :: this
    real(kind=prec), intent(in) :: s(1:this%nvar)
    real(kind=prec), intent(in) :: nhat(1:2)

    Return Value real(kind=prec)(1:this%nvar)

procedure, public :: hbc3d_NoNormalFlow => hbc3d_Generic_Model

  • public pure function hbc3d_Generic_Model(this, s, nhat) result(exts)

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(in) :: this
    real(kind=prec), intent(in) :: s(1:this%nvar)
    real(kind=prec), intent(in) :: nhat(1:3)

    Return Value real(kind=prec)(1:this%nvar)

procedure, public :: hbc3d_Prescribed => hbc3d_Prescribed_Model

  • public pure function hbc3d_Prescribed_Model(this, x, t) result(exts)

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(in) :: this
    real(kind=prec), intent(in) :: x(1:3)
    real(kind=prec), intent(in) :: t

    Return Value real(kind=prec)(1:this%nvar)

procedure, public :: hbc3d_Radiation => hbc3d_Generic_Model

  • public pure function hbc3d_Generic_Model(this, s, nhat) result(exts)

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(in) :: this
    real(kind=prec), intent(in) :: s(1:this%nvar)
    real(kind=prec), intent(in) :: nhat(1:3)

    Return Value real(kind=prec)(1:this%nvar)

procedure, public :: pbc1d_NoNormalFlow => pbc1d_Generic_Model

  • public pure function pbc1d_Generic_Model(this, dsdx, nhat) result(extDsdx)

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(in) :: this
    real(kind=prec), intent(in) :: dsdx(1:this%nvar)
    real(kind=prec), intent(in) :: nhat

    Return Value real(kind=prec)(1:this%nvar)

procedure, public :: pbc1d_Prescribed => pbc1d_Prescribed_Model

  • public pure function pbc1d_Prescribed_Model(this, x, t) result(extDsdx)

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(in) :: this
    real(kind=prec), intent(in) :: x
    real(kind=prec), intent(in) :: t

    Return Value real(kind=prec)(1:this%nvar)

procedure, public :: pbc1d_Radiation => pbc1d_Generic_Model

  • public pure function pbc1d_Generic_Model(this, dsdx, nhat) result(extDsdx)

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(in) :: this
    real(kind=prec), intent(in) :: dsdx(1:this%nvar)
    real(kind=prec), intent(in) :: nhat

    Return Value real(kind=prec)(1:this%nvar)

procedure, public :: pbc2d_NoNormalFlow => pbc2d_Generic_Model

  • public pure function pbc2d_Generic_Model(this, dsdx, nhat) result(extDsdx)

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(in) :: this
    real(kind=prec), intent(in) :: dsdx(1:this%nvar,1:2)
    real(kind=prec), intent(in) :: nhat(1:2)

    Return Value real(kind=prec)(1:this%nvar,1:2)

procedure, public :: pbc2d_Prescribed => pbc2d_Prescribed_Model

  • public pure function pbc2d_Prescribed_Model(this, x, t) result(extDsdx)

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(in) :: this
    real(kind=prec), intent(in) :: x(1:2)
    real(kind=prec), intent(in) :: t

    Return Value real(kind=prec)(1:this%nvar,1:2)

procedure, public :: pbc2d_Radiation => pbc2d_Generic_Model

  • public pure function pbc2d_Generic_Model(this, dsdx, nhat) result(extDsdx)

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(in) :: this
    real(kind=prec), intent(in) :: dsdx(1:this%nvar,1:2)
    real(kind=prec), intent(in) :: nhat(1:2)

    Return Value real(kind=prec)(1:this%nvar,1:2)

procedure, public :: pbc3d_NoNormalFlow => pbc3d_Generic_Model

  • public pure function pbc3d_Generic_Model(this, dsdx, nhat) result(extDsdx)

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(in) :: this
    real(kind=prec), intent(in) :: dsdx(1:this%nvar,1:3)
    real(kind=prec), intent(in) :: nhat(1:3)

    Return Value real(kind=prec)(1:this%nvar,1:3)

procedure, public :: pbc3d_Prescribed => pbc3d_Prescribed_Model

  • public pure function pbc3d_Prescribed_Model(this, x, t) result(extDsdx)

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(in) :: this
    real(kind=prec), intent(in) :: x(1:3)
    real(kind=prec), intent(in) :: t

    Return Value real(kind=prec)(1:this%nvar,1:3)

procedure, public :: pbc3d_Radiation => pbc3d_Generic_Model

  • public pure function pbc3d_Generic_Model(this, dsdx, nhat) result(extDsdx)

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(in) :: this
    real(kind=prec), intent(in) :: dsdx(1:this%nvar,1:3)
    real(kind=prec), intent(in) :: nhat(1:3)

    Return Value real(kind=prec)(1:this%nvar,1:3)

procedure, public :: riemannflux1d => riemannflux1d_Model

  • public pure function riemannflux1d_Model(this, sL, sR, dsdx, nhat) result(flux)

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(in) :: this
    real(kind=prec), intent(in) :: sL(1:this%nvar)
    real(kind=prec), intent(in) :: sR(1:this%nvar)
    real(kind=prec), intent(in) :: dsdx(1:this%nvar)
    real(kind=prec), intent(in) :: nhat

    Return Value real(kind=prec)(1:this%nvar)

procedure, public :: riemannflux2d => riemannflux2d_Model

  • public pure function riemannflux2d_Model(this, sL, sR, dsdx, nhat) result(flux)

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(in) :: this
    real(kind=prec), intent(in) :: sL(1:this%nvar)
    real(kind=prec), intent(in) :: sR(1:this%nvar)
    real(kind=prec), intent(in) :: dsdx(1:this%nvar,1:2)
    real(kind=prec), intent(in) :: nhat(1:2)

    Return Value real(kind=prec)(1:this%nvar)

procedure, public :: riemannflux3d => riemannflux3d_Model

  • public pure function riemannflux3d_Model(this, sL, sR, dsdx, nhat) result(flux)

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(in) :: this
    real(kind=prec), intent(in) :: sL(1:this%nvar)
    real(kind=prec), intent(in) :: sR(1:this%nvar)
    real(kind=prec), intent(in) :: dsdx(1:this%nvar,1:3)
    real(kind=prec), intent(in) :: nhat(1:3)

    Return Value real(kind=prec)(1:this%nvar)

procedure, public :: source1d => source1d_Model

  • public pure function source1d_Model(this, s, dsdx) result(source)

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(in) :: this
    real(kind=prec), intent(in) :: s(1:this%nvar)
    real(kind=prec), intent(in) :: dsdx(1:this%nvar)

    Return Value real(kind=prec)(1:this%nvar)

procedure, public :: source2d => source2d_Model

  • public pure function source2d_Model(this, s, dsdx) result(source)

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(in) :: this
    real(kind=prec), intent(in) :: s(1:this%nvar)
    real(kind=prec), intent(in) :: dsdx(1:this%nvar,1:2)

    Return Value real(kind=prec)(1:this%nvar)

procedure, public :: source3d => source3d_Model

  • public pure function source3d_Model(this, s, dsdx) result(source)

    Arguments

    TypeIntentOptionalAttributesName
    class(Model), intent(in) :: this
    real(kind=prec), intent(in) :: s(1:this%nvar)
    real(kind=prec), intent(in) :: dsdx(1:this%nvar,1:3)

    Return Value real(kind=prec)(1:this%nvar)

Source Code

  type,abstract :: Model

    ! Time integration attributes
    procedure(SELF_timeIntegrator),pointer :: timeIntegrator => Euler_timeIntegrator
    real(prec) :: dt
    real(prec) :: t
    integer :: ioIterate = 0
    logical :: gradient_enabled = .false.
    logical :: prescribed_bcs_enabled = .true.
    logical :: tecplot_enabled = .true.
    integer :: nvar
    ! Standard Diagnostics
    real(prec) :: entropy ! Mathematical entropy function for the model

  contains

    procedure :: IncrementIOCounter

    procedure :: PrintType => PrintType_Model

    procedure :: SetNumberOfVariables => SetNumberOfVariables_Model

    procedure :: AdditionalInit => AdditionalInit_Model
    procedure :: AdditionalFree => AdditionalFree_Model
    procedure :: AdditionalOutput => AdditionalOutput_Model

    procedure :: ForwardStep => ForwardStep_Model

    procedure :: Euler_timeIntegrator

    ! Runge-Kutta methods
    procedure :: LowStorageRK2_timeIntegrator
    procedure(UpdateGRK),deferred :: UpdateGRK2

    procedure :: LowStorageRK3_timeIntegrator
    procedure(UpdateGRK),deferred :: UpdateGRK3

    procedure :: LowStorageRK4_timeIntegrator
    procedure(UpdateGRK),deferred :: UpdateGRK4

    procedure :: PreTendency => PreTendency_Model
    procedure :: entropy_func => entropy_func_Model

    procedure :: flux1D => flux1d_Model
    procedure :: flux2D => flux2d_Model
    procedure :: flux3D => flux3d_Model

    procedure :: riemannflux1d => riemannflux1d_Model
    procedure :: riemannflux2d => riemannflux2d_Model
    procedure :: riemannflux3d => riemannflux3d_Model

    procedure :: source1d => source1d_Model
    procedure :: source2d => source2d_Model
    procedure :: source3d => source3d_Model

    ! Boundary condition functions (hyperbolic)
    procedure :: hbc1d_Prescribed => hbc1d_Prescribed_Model
    procedure :: hbc1d_Radiation => hbc1d_Generic_Model
    procedure :: hbc1d_NoNormalFlow => hbc1d_Generic_Model
    procedure :: hbc2d_Prescribed => hbc2d_Prescribed_Model
    procedure :: hbc2d_Radiation => hbc2d_Generic_Model
    procedure :: hbc2d_NoNormalFlow => hbc2d_Generic_Model
    procedure :: hbc3d_Prescribed => hbc3d_Prescribed_Model
    procedure :: hbc3d_Radiation => hbc3d_Generic_Model
    procedure :: hbc3d_NoNormalFlow => hbc3d_Generic_Model

    ! Boundary condition functions (parabolic)
    procedure :: pbc1d_Prescribed => pbc1d_Prescribed_Model
    procedure :: pbc1d_Radiation => pbc1d_Generic_Model
    procedure :: pbc1d_NoNormalFlow => pbc1d_Generic_Model
    procedure :: pbc2d_Prescribed => pbc2d_Prescribed_Model
    procedure :: pbc2d_Radiation => pbc2d_Generic_Model
    procedure :: pbc2d_NoNormalFlow => pbc2d_Generic_Model
    procedure :: pbc3d_Prescribed => pbc3d_Prescribed_Model
    procedure :: pbc3d_Radiation => pbc3d_Generic_Model
    procedure :: pbc3d_NoNormalFlow => pbc3d_Generic_Model

    procedure :: ReportEntropy => ReportEntropy_Model
    procedure :: ReportMetrics => ReportMetrics_Model
    procedure :: ReportUserMetrics => ReportUserMetrics_Model
    procedure :: CalculateEntropy => CalculateEntropy_Model

    procedure(UpdateSolution),deferred :: UpdateSolution
    procedure(CalculateTendency),deferred :: CalculateTendency
    procedure(ReadModel),deferred :: ReadModel
    procedure(WriteModel),deferred :: WriteModel
    procedure(WriteTecplot),deferred :: WriteTecplot

    generic :: SetTimeIntegrator => SetTimeIntegrator_withChar
    procedure,private :: SetTimeIntegrator_withChar

    procedure :: SetSimulationTime
    procedure :: GetSimulationTime

  endtype Model