integrand_lcce Derived Type

type, public, extends(integrand_tester_object) :: integrand_lcce

type~~integrand_lcce~~InheritsGraph type~integrand_lcce integrand_lcce type~integrand_tester_object integrand_tester_object type~integrand_tester_object->type~integrand_lcce type~integrand_object integrand_object type~integrand_object->type~integrand_tester_object
Help


The linear constant coefficient equation field.

It is a FOODIE integrand class concrete extension.

ODE

The linear constant coefficient equation is a linear pure ODE and it can be written as:

$$\begin{matrix} U_t = R(U) \\ U = \begin{bmatrix} y \end{bmatrix}\;\;\; R(U) = \begin{bmatrix} a*y + b \end{bmatrix} \end{matrix}$$

The coefficent a,b are constant with \(a \ne 0\). The exact solution is

$$ U(t) = (U_0 + \frac{b}{a}) e ^{a(t-t_0)} - \frac{b){a} $$

where \(t_0\) is the initial time of integration.

Bibliography

[1] ORDINARY DIFFERENTIAL EQUATIONS, Gabriel Nagy, 2017.

Inherited By

type~~integrand_lcce~~InheritedByGraph type~integrand_lcce integrand_lcce type~test_object test_object type~integrand_lcce->type~test_object lcce_0
Help

Source Code


Components

TypeVisibility AttributesNameInitial
real(kind=R_P), private :: U =0._R_P

Integrand (state) variable.

real(kind=R_P), private :: U0 =0._R_P

Integrand initial state.

real(kind=R_P), private :: a =0._R_P

a constant.

real(kind=R_P), private :: b =0._R_P

b constant.


Type-Bound Procedures

generic, public :: add_fast => integrand_add_integrand_fast

Overloading add_fast method.

procedure, public, pass(lhs) :: assign_integrand

= operator.

procedure, public, pass(lhs) :: assign_real

= real operator.

  • private pure subroutine assign_real(lhs, rhs)

    = real operator.

    Arguments

    Type IntentOptional AttributesName
    class(integrand_lcce), intent(inout) :: lhs

    Left hand side.

    real(kind=R_P), intent(in) :: rhs(1:)

    Right hand side.

generic, public :: assignment(=) => assign_integrand, assign_real

Overloading = assignament.

  • private pure subroutine assign_integrand(lhs, rhs)

    = operator.

    Arguments

    Type IntentOptional AttributesName
    class(integrand_oscillation), intent(inout) :: lhs

    Left hand side.

    class(integrand_object), intent(in) :: rhs

    Right hand side.

  • private pure subroutine assign_real(lhs, rhs)

    = real operator.

    Arguments

    Type IntentOptional AttributesName
    class(integrand_oscillation), intent(inout) :: lhs

    Left hand side.

    real(kind=R_P), intent(in) :: rhs(1:)

    Right hand side.

procedure, public, pass(self) :: description

Return an informative description of the test.

  • private pure function description(self, prefix) result(desc)

    Return informative integrator description.

    Arguments

    Type IntentOptional AttributesName
    class(integrand_lcce), intent(in) :: self

    Integrand.

    character(len=*), intent(in), optional :: prefix

    Prefixing string.

    Return Value character(len=:), allocatable

    Description.

procedure, public, pass(self) :: error

Return error.

  • private pure function error(self, t, t0, U0)

    Return error.

    Arguments

    Type IntentOptional AttributesName
    class(integrand_lcce), intent(in) :: self

    Integrand.

    real(kind=R_P), intent(in) :: t

    Time.

    real(kind=R_P), intent(in), optional :: t0

    Initial time.

    class(integrand_object), intent(in), optional :: U0

    Initial conditions.

    Return Value real(kind=R_P), allocatable, (:)

    Error.

procedure, public, pass(self) :: exact_solution

Return exact solution.

  • private pure function exact_solution(self, t, t0, U0) result(exact)

    Return exact solution.

    Arguments

    Type IntentOptional AttributesName
    class(integrand_lcce), intent(in) :: self

    Integrand.

    real(kind=R_P), intent(in) :: t

    Time.

    real(kind=R_P), intent(in), optional :: t0

    Initial time.

    class(integrand_object), intent(in), optional :: U0

    Initial conditions.

    Return Value real(kind=R_P), allocatable, (:)

    Exact solution.

procedure, public, pass(self) :: export_tecplot

Export integrand to Tecplot file.

  • private subroutine export_tecplot(self, file_name, t, scheme, close_file, with_exact_solution, U0)

    Export integrand to Tecplot file.

    Arguments

    Type IntentOptional AttributesName
    class(integrand_lcce), intent(in) :: self

    Advection field.

    character(len=*), intent(in), optional :: file_name

    File name.

    real(kind=R_P), intent(in), optional :: t

    Time.

    character(len=*), intent(in), optional :: scheme

    Scheme used to integrate integrand.

    logical, intent(in), optional :: close_file

    Flag for closing file.

    logical, intent(in), optional :: with_exact_solution

    Flag for export also exact solution.

    class(integrand_object), intent(in), optional :: U0

    Initial conditions.

procedure, public, pass(self) :: initialize

Initialize field.

  • private pure subroutine initialize(self, Dt)

    Initialize integrand.

    Arguments

    Type IntentOptional AttributesName
    class(integrand_lcce), intent(inout) :: self

    Integrand.

    real(kind=R_P), intent(in) :: Dt

    Time step.

procedure, public, pass(lhs) :: integrand_add_integrand

+ operator.

  • private pure function integrand_add_integrand(lhs, rhs) result(opr)

    + operator.

    Arguments

    Type IntentOptional AttributesName
    class(integrand_lcce), intent(in) :: lhs

    Left hand side.

    class(integrand_object), intent(in) :: rhs

    Right hand side.

    Return Value real(kind=R_P), allocatable, (:)

    Operator result.

procedure, public, pass(opr) :: integrand_add_integrand_fast

+ fast operator.

procedure, public, pass(lhs) :: integrand_add_real

+ real operator.

  • private pure function integrand_add_real(lhs, rhs) result(opr)

    + real operator.

    Arguments

    Type IntentOptional AttributesName
    class(integrand_lcce), intent(in) :: lhs

    Left hand side.

    real(kind=R_P), intent(in) :: rhs(1:)

    Right hand side.

    Return Value real(kind=R_P), allocatable, (:)

    Operator result.

procedure, public, pass(self) :: integrand_dimension

Return integrand dimension.

  • private pure function integrand_dimension(self)

    return integrand dimension.

    Arguments

    Type IntentOptional AttributesName
    class(integrand_lcce), intent(in) :: self

    integrand.

    Return Value integer(kind=I_P)

    integrand dimension.

procedure, public, pass(lhs) :: integrand_multiply_integrand

* operator.

  • private pure function integrand_multiply_integrand(lhs, rhs) result(opr)

    * operator.

    Arguments

    Type IntentOptional AttributesName
    class(integrand_lcce), intent(in) :: lhs

    Left hand side.

    class(integrand_object), intent(in) :: rhs

    Right hand side.

    Return Value real(kind=R_P), allocatable, (:)

    Operator result.

procedure, public, pass(opr) :: integrand_multiply_integrand_fast

* fast operator.

procedure, public, pass(lhs) :: integrand_multiply_real

* real operator.

  • private pure function integrand_multiply_real(lhs, rhs) result(opr)

    * real_scalar operator.

    Arguments

    Type IntentOptional AttributesName
    class(integrand_lcce), intent(in) :: lhs

    Left hand side.

    real(kind=R_P), intent(in) :: rhs(1:)

    Right hand side.

    Return Value real(kind=R_P), allocatable, (:)

    Operator result.

procedure, public, pass(lhs) :: integrand_multiply_real_scalar

* real_scalar operator.

  • private pure function integrand_multiply_real_scalar(lhs, rhs) result(opr)

    * real_scalar operator.

    Arguments

    Type IntentOptional AttributesName
    class(integrand_lcce), intent(in) :: lhs

    Left hand side.

    real(kind=R_P), intent(in) :: rhs

    Right hand side.

    Return Value real(kind=R_P), allocatable, (:)

    Operator result.

procedure, public, pass(opr) :: integrand_multiply_real_scalar_fast

* real_scalar fast operator.

procedure, public, pass(lhs) :: integrand_sub_integrand

- operator.

  • private pure function integrand_sub_integrand(lhs, rhs) result(opr)

    - operator.

    Arguments

    Type IntentOptional AttributesName
    class(integrand_lcce), intent(in) :: lhs

    Left hand side.

    class(integrand_object), intent(in) :: rhs

    Right hand side.

    Return Value real(kind=R_P), allocatable, (:)

    Operator result.

procedure, public, pass(lhs) :: integrand_sub_real

- real operator.

  • private pure function integrand_sub_real(lhs, rhs) result(opr)

    - real operator.

    Arguments

    Type IntentOptional AttributesName
    class(integrand_lcce), intent(in) :: lhs

    Left hand side.

    real(kind=R_P), intent(in) :: rhs(1:)

    Right hand side.

    Return Value real(kind=R_P), allocatable, (:)

    Operator result.

procedure, public, pass(opr) :: integrand_subtract_integrand_fast

- fast operator.

procedure, public, pass(lhs) :: local_error

||integrand_lcce - integrand_lcce|| operator.

  • private pure function local_error(lhs, rhs) result(error)

    Estimate local truncation error between 2 oscillation approximations.

    Arguments

    Type IntentOptional AttributesName
    class(integrand_lcce), intent(in) :: lhs

    Left hand side.

    class(integrand_object), intent(in) :: rhs

    Right hand side.

    Return Value real(kind=R_P)

    Error estimation.

Overloading multiply_fast method.

Overloading * operator.

  • private pure function integrand_multiply_integrand(lhs, rhs) result(opr)

    * operator.

    Arguments

    Type IntentOptional AttributesName
    class(integrand_oscillation), intent(in) :: lhs

    Left hand side.

    class(integrand_object), intent(in) :: rhs

    Right hand side.

    Return Value real(kind=R_P), allocatable, (:)

    Operator result.

  • private pure function integrand_multiply_real(lhs, rhs) result(opr)

    * real_scalar operator.

    Arguments

    Type IntentOptional AttributesName
    class(integrand_oscillation), intent(in) :: lhs

    Left hand side.

    real(kind=R_P), intent(in) :: rhs(1:)

    Right hand side.

    Return Value real(kind=R_P), allocatable, (:)

    Operator result.

  • private pure function real_multiply_integrand(lhs, rhs) result(opr)

    real_scalar * operator.

    Arguments

    Type IntentOptional AttributesName
    real(kind=R_P), intent(in) :: lhs(1:)

    Left hand side.

    class(integrand_oscillation), intent(in) :: rhs

    Right hand side.

    Return Value real(kind=R_P), allocatable, (:)

    Operator result.

  • private pure function integrand_multiply_real_scalar(lhs, rhs) result(opr)

    * real_scalar operator.

    Arguments

    Type IntentOptional AttributesName
    class(integrand_oscillation), intent(in) :: lhs

    Left hand side.

    real(kind=R_P), intent(in) :: rhs

    Right hand side.

    Return Value real(kind=R_P), allocatable, (:)

    Operator result.

  • private pure function real_scalar_multiply_integrand(lhs, rhs) result(opr)

    real_scalar * operator.

    Arguments

    Type IntentOptional AttributesName
    real(kind=R_P), intent(in) :: lhs

    Left hand side.

    class(integrand_oscillation), intent(in) :: rhs

    Right hand side.

    Return Value real(kind=R_P), allocatable, (:)

    Operator result.

Overloading + operator.

  • private pure function integrand_add_integrand(lhs, rhs) result(opr)

    + operator.

    Arguments

    Type IntentOptional AttributesName
    class(integrand_oscillation), intent(in) :: lhs

    Left hand side.

    class(integrand_object), intent(in) :: rhs

    Right hand side.

    Return Value real(kind=R_P), allocatable, (:)

    Operator result.

  • private pure function integrand_add_real(lhs, rhs) result(opr)

    + real operator.

    Arguments

    Type IntentOptional AttributesName
    class(integrand_oscillation), intent(in) :: lhs

    Left hand side.

    real(kind=R_P), intent(in) :: rhs(1:)

    Right hand side.

    Return Value real(kind=R_P), allocatable, (:)

    Operator result.

  • private pure function real_add_integrand(lhs, rhs) result(opr)

    real + operator.

    Arguments

    Type IntentOptional AttributesName
    real(kind=R_P), intent(in) :: lhs(1:)

    Left hand side.

    class(integrand_oscillation), intent(in) :: rhs

    Left hand side.

    Return Value real(kind=R_P), allocatable, (:)

    Operator result.

Overloading - operator.

  • private pure function integrand_sub_integrand(lhs, rhs) result(opr)

    - operator.

    Arguments

    Type IntentOptional AttributesName
    class(integrand_oscillation), intent(in) :: lhs

    Left hand side.

    class(integrand_object), intent(in) :: rhs

    Right hand side.

    Return Value real(kind=R_P), allocatable, (:)

    Operator result.

  • private pure function integrand_sub_real(lhs, rhs) result(opr)

    - real operator.

    Arguments

    Type IntentOptional AttributesName
    class(integrand_oscillation), intent(in) :: lhs

    Left hand side.

    real(kind=R_P), intent(in) :: rhs(1:)

    Right hand side.

    Return Value real(kind=R_P), allocatable, (:)

    Operator result.

  • private pure function real_sub_integrand(lhs, rhs) result(opr)

    real - operator.

    Arguments

    Type IntentOptional AttributesName
    real(kind=R_P), intent(in) :: lhs(1:)

    Left hand side.

    class(integrand_oscillation), intent(in) :: rhs

    Left hand side.

    Return Value real(kind=R_P), allocatable, (:)

    Operator result.

generic, public :: operator(.lterror.) => local_error

Estimate local truncation error.

  • private pure function local_error(lhs, rhs) result(error)

    Estimate local truncation error between 2 oscillation approximations.

    Arguments

    Type IntentOptional AttributesName
    class(integrand_oscillation), intent(in) :: lhs

    Left hand side.

    class(integrand_object), intent(in) :: rhs

    Right hand side.

    Return Value real(kind=R_P)

    Error estimation.

procedure, public, pass(self) :: output

Extract integrand state field.

  • private pure function output(self) result(state)

    Extract integrand state field.

    Arguments

    Type IntentOptional AttributesName
    class(integrand_lcce), intent(in) :: self

    Integrand.

    Return Value real(kind=R_P)

    State.

procedure, public, pass(self) :: parse_cli

Initialize from command line interface.

  • private subroutine parse_cli(self, cli)

    Initialize from command line interface.

    Arguments

    Type IntentOptional AttributesName
    class(integrand_lcce), intent(inout) :: self

    Advection field.

    type(command_line_interface), intent(inout) :: cli

    Command line interface handler.

procedure, public, pass(rhs) :: real_add_integrand

real + operator.

  • private pure function real_add_integrand(lhs, rhs) result(opr)

    real + operator.

    Arguments

    Type IntentOptional AttributesName
    real(kind=R_P), intent(in) :: lhs(1:)

    Left hand side.

    class(integrand_lcce), intent(in) :: rhs

    Left hand side.

    Return Value real(kind=R_P), allocatable, (:)

    Operator result.

procedure, public, pass(rhs) :: real_multiply_integrand

real * operator.

  • private pure function real_multiply_integrand(lhs, rhs) result(opr)

    real_scalar * operator.

    Arguments

    Type IntentOptional AttributesName
    real(kind=R_P), intent(in) :: lhs(1:)

    Left hand side.

    class(integrand_lcce), intent(in) :: rhs

    Right hand side.

    Return Value real(kind=R_P), allocatable, (:)

    Operator result.

procedure, public, pass(rhs) :: real_scalar_multiply_integrand

real_scalar * operator.

  • private pure function real_scalar_multiply_integrand(lhs, rhs) result(opr)

    real_scalar * operator.

    Arguments

    Type IntentOptional AttributesName
    real(kind=R_P), intent(in) :: lhs

    Left hand side.

    class(integrand_lcce), intent(in) :: rhs

    Right hand side.

    Return Value real(kind=R_P), allocatable, (:)

    Operator result.

procedure, public, pass(rhs) :: real_sub_integrand

real - operator.

  • private pure function real_sub_integrand(lhs, rhs) result(opr)

    real - operator.

    Arguments

    Type IntentOptional AttributesName
    real(kind=R_P), intent(in) :: lhs(1:)

    Left hand side.

    class(integrand_lcce), intent(in) :: rhs

    Left hand side.

    Return Value real(kind=R_P), allocatable, (:)

    Operator result.

procedure, public, nopass :: set_cli

Set command line interface.

  • private subroutine set_cli(cli)

    Set command line interface.

    Arguments

    Type IntentOptional AttributesName
    type(command_line_interface), intent(inout) :: cli

    Command line interface handler.

generic, public :: subtract_fast => integrand_subtract_integrand_fast

Overloading subtract_fast method.

procedure, public, pass(self) :: t => dU_dt

Time derivative, residuals.

  • private pure function dU_dt(self, t) result(dState_dt)

    Time derivative of field.

    Arguments

    Type IntentOptional AttributesName
    class(integrand_lcce), intent(in) :: self

    Integrand.

    real(kind=R_P), intent(in), optional :: t

    Time.

    Return Value real(kind=R_P), allocatable, (:)

    Integrand time derivative.

procedure, public, pass(self) :: t_fast

Time derivative, residuals, fast mode.

  • private subroutine t_fast(self, t)

    Time derivative function of integrand class, i.e. the residuals function. Fast mode acting directly on self.

    Arguments

    Type IntentOptional AttributesName
    class(integrand_lcce), intent(inout) :: self

    Integrand.

    real(kind=R_P), intent(in), optional :: t

    Time.

Source Code

type, extends(integrand_tester_object) :: integrand_lcce
   !< The linear constant coefficient equation field.
   !<
   !< It is a FOODIE integrand class concrete extension.
   !<
   !<### ODE
   !<The linear constant coefficient equation is a linear pure ODE and it can be written as:
   !<
   !<$$\begin{matrix}
   !< U_t = R(U)  \\
   !< U = \begin{bmatrix}
   !< y
   !< \end{bmatrix}\;\;\;
   !< R(U) = \begin{bmatrix}
   !< a*y + b
   !< \end{bmatrix}
   !<\end{matrix}$$
   !<
   !<The coefficent *a,b* are constant with \(a \ne 0\). The exact solution is
   !<
   !<$$ U(t) = (U_0 + \frac{b}{a}) e ^{a(t-t_0)} - \frac{b){a} $$
   !<
   !< where \(t_0\) is the initial time of integration.
   !<
   !<#### Bibliography
   !<
   !<[1] *ORDINARY DIFFERENTIAL EQUATIONS*, Gabriel Nagy, 2017.
   !<
   private
   real(R_P) :: a=0._R_P  !< *a* constant.
   real(R_P) :: b=0._R_P  !< *b* constant.
   real(R_P) :: U=0._R_P  !< Integrand (state) variable.
   real(R_P) :: U0=0._R_P !< Integrand initial state.
   contains
      ! auxiliary methods
      procedure, pass(self), public :: output         !< Extract integrand state field.
      ! integrand_tester_object deferred methods
      procedure, pass(self), public :: description    !< Return an informative description of the test.
      procedure, pass(self), public :: error          !< Return error.
      procedure, pass(self), public :: exact_solution !< Return exact solution.
      procedure, pass(self), public :: export_tecplot !< Export integrand to Tecplot file.
      procedure, pass(self), public :: initialize     !< Initialize field.
      procedure, pass(self), public :: parse_cli      !< Initialize from command line interface.
      procedure, nopass,     public :: set_cli        !< Set command line interface.
      ! integrand_object deferred methods
      procedure, pass(self), public :: integrand_dimension !< Return integrand dimension.
      procedure, pass(self), public :: t => dU_dt          !< Time derivative, residuals.
      ! operators
      procedure, pass(lhs), public :: local_error !<`||integrand_lcce - integrand_lcce||` operator.
      ! +
      procedure, pass(lhs), public :: integrand_add_integrand !< `+` operator.
      procedure, pass(lhs), public :: integrand_add_real      !< `+ real` operator.
      procedure, pass(rhs), public :: real_add_integrand      !< `real +` operator.
      ! *
      procedure, pass(lhs), public :: integrand_multiply_integrand   !< `*` operator.
      procedure, pass(lhs), public :: integrand_multiply_real        !< `* real` operator.
      procedure, pass(rhs), public :: real_multiply_integrand        !< `real *` operator.
      procedure, pass(lhs), public :: integrand_multiply_real_scalar !< `* real_scalar` operator.
      procedure, pass(rhs), public :: real_scalar_multiply_integrand !< `real_scalar *` operator.
      ! -
      procedure, pass(lhs), public :: integrand_sub_integrand !< `-` operator.
      procedure, pass(lhs), public :: integrand_sub_real      !< `- real` operator.
      procedure, pass(rhs), public :: real_sub_integrand      !< `real -` operator.
      ! =
      procedure, pass(lhs), public :: assign_integrand !< `=` operator.
      procedure, pass(lhs), public :: assign_real      !< `= real` operator.
      ! override fast operators
      procedure, pass(self), public :: t_fast                              !< Time derivative, residuals, fast mode.
      procedure, pass(opr),  public :: integrand_add_integrand_fast        !< `+` fast operator.
      procedure, pass(opr),  public :: integrand_multiply_integrand_fast   !< `*` fast operator.
      procedure, pass(opr),  public :: integrand_multiply_real_scalar_fast !< `* real_scalar` fast operator.
      procedure, pass(opr),  public :: integrand_subtract_integrand_fast   !< `-` fast operator.
endtype integrand_lcce