pf_level_t Derived Type

type, public :: pf_level_t

Data type of a PFASST level


Inherits

type~~pf_level_t~~InheritsGraph type~pf_level_t pf_level_t type~pf_user_level_t pf_user_level_t type~pf_level_t->type~pf_user_level_t ulevel type~pf_encap_t pf_encap_t type~pf_level_t->type~pf_encap_t Q, pQ, R, I, Fflt, tauQ, pFflt, q0, qend, F, pF type~pf_sdcmats_t pf_sdcmats_t type~pf_level_t->type~pf_sdcmats_t sdcmats type~pf_sweeper_t pf_sweeper_t type~pf_user_level_t->type~pf_sweeper_t sweeper type~pf_factory_t pf_factory_t type~pf_user_level_t->type~pf_factory_t factory type~pf_stepper_t pf_stepper_t type~pf_user_level_t->type~pf_stepper_t stepper

Inherited by

type~~pf_level_t~~InheritedByGraph type~pf_level_t pf_level_t type~pf_pfasst_t pf_pfasst_t type~pf_pfasst_t->type~pf_level_t levels

Contents

Source Code


Components

TypeVisibility AttributesNameInitial
integer, public :: index =-1

level number (1 is the coarsest)

integer, public :: nnodes =-1

number of sdc nodes

integer, public :: nsteps_rk =-1

number of rk steps to perform

integer, public :: nsweeps =-1

number of sdc sweeps to perform

integer, public :: nsweeps_pred =-1

number of coarse sdc sweeps to perform predictor in predictor

logical, public :: Finterp =.false.

interpolate functions instead of solutions

integer, public :: mpibuflen =-1

size of solution in pfdp units

real(kind=pfdp), public :: error

holds the user defined error

real(kind=pfdp), public :: residual

holds the user defined residual

real(kind=pfdp), public :: residual_rel

holds the user defined relative residual (scaled by solution magnitude)

class(pf_user_level_t), public, allocatable:: ulevel

user customized level info

real(kind=pfdp), public, allocatable:: send(:)

Simple data storage at each level send buffer recv buffer list of SDC nodes time restriction matrix time interpolation matrix

real(kind=pfdp), public, allocatable:: recv(:)

Simple data storage at each level send buffer recv buffer list of SDC nodes time restriction matrix time interpolation matrix

real(kind=pfdp), public, allocatable:: nodes(:)

Simple data storage at each level send buffer recv buffer list of SDC nodes time restriction matrix time interpolation matrix

real(kind=pfdp), public, allocatable:: rmat(:,:)

Simple data storage at each level send buffer recv buffer list of SDC nodes time restriction matrix time interpolation matrix

real(kind=pfdp), public, allocatable:: tmat(:,:)

Simple data storage at each level send buffer recv buffer list of SDC nodes time restriction matrix time interpolation matrix

integer, public, allocatable:: nflags(:)

sdc node flags

class(pf_encap_t), public, allocatable:: Q(:)

Solution variable storage solution at sdc nodes unknowns at sdc nodes, previous sweep full residuals 0 to node integrals functions values at sdc nodes (flat) fas correction in Q form functions at sdc nodes, previous sweep (flat) initial condition solution at end time

class(pf_encap_t), public, allocatable:: pQ(:)

Solution variable storage solution at sdc nodes unknowns at sdc nodes, previous sweep full residuals 0 to node integrals functions values at sdc nodes (flat) fas correction in Q form functions at sdc nodes, previous sweep (flat) initial condition solution at end time

class(pf_encap_t), public, allocatable:: R(:)

Solution variable storage solution at sdc nodes unknowns at sdc nodes, previous sweep full residuals 0 to node integrals functions values at sdc nodes (flat) fas correction in Q form functions at sdc nodes, previous sweep (flat) initial condition solution at end time

class(pf_encap_t), public, allocatable:: I(:)

Solution variable storage solution at sdc nodes unknowns at sdc nodes, previous sweep full residuals 0 to node integrals functions values at sdc nodes (flat) fas correction in Q form functions at sdc nodes, previous sweep (flat) initial condition solution at end time

class(pf_encap_t), public, allocatable:: Fflt(:)

Solution variable storage solution at sdc nodes unknowns at sdc nodes, previous sweep full residuals 0 to node integrals functions values at sdc nodes (flat) fas correction in Q form functions at sdc nodes, previous sweep (flat) initial condition solution at end time

class(pf_encap_t), public, allocatable:: tauQ(:)

Solution variable storage solution at sdc nodes unknowns at sdc nodes, previous sweep full residuals 0 to node integrals functions values at sdc nodes (flat) fas correction in Q form functions at sdc nodes, previous sweep (flat) initial condition solution at end time

class(pf_encap_t), public, allocatable:: pFflt(:)

Solution variable storage solution at sdc nodes unknowns at sdc nodes, previous sweep full residuals 0 to node integrals functions values at sdc nodes (flat) fas correction in Q form functions at sdc nodes, previous sweep (flat) initial condition solution at end time

class(pf_encap_t), public, allocatable:: q0

Solution variable storage solution at sdc nodes unknowns at sdc nodes, previous sweep full residuals 0 to node integrals functions values at sdc nodes (flat) fas correction in Q form functions at sdc nodes, previous sweep (flat) initial condition solution at end time

class(pf_encap_t), public, allocatable:: qend

Solution variable storage solution at sdc nodes unknowns at sdc nodes, previous sweep full residuals 0 to node integrals functions values at sdc nodes (flat) fas correction in Q form functions at sdc nodes, previous sweep (flat) initial condition solution at end time

class(pf_encap_t), public, pointer:: F(:,:)

Function storage functions values at sdc nodes functions at sdc nodes, previous sweep

class(pf_encap_t), public, pointer:: pF(:,:)

Function storage functions values at sdc nodes functions at sdc nodes, previous sweep

integer, public, allocatable:: shape(:)

user defined shape array

type(pf_sdcmats_t), public, allocatable:: sdcmats
logical, public :: allocated =.false.

Source Code

  type :: pf_level_t
     !  level parameters set by the pfasst_t values
     integer  :: index        = -1   !! level number (1 is the coarsest)
     integer  :: nnodes       = -1   !! number of sdc nodes
     integer  :: nsteps_rk    = -1   !! number of rk steps to perform
     integer  :: nsweeps      = -1   !! number of sdc sweeps to perform
     integer  :: nsweeps_pred = -1      !! number of coarse sdc sweeps to perform predictor in predictor
     logical     :: Finterp = .false.   !! interpolate functions instead of solutions

     !  Mandatory level parameter
     integer  :: mpibuflen    = -1   !! size of solution in pfdp units



     !  Diagnostics
     real(pfdp)  :: error            !! holds the user defined error
     real(pfdp)  :: residual         !! holds the user defined residual
     real(pfdp)  :: residual_rel     !! holds the user defined relative residual (scaled by solution magnitude)

     class(pf_user_level_t), allocatable :: ulevel  !!  user customized level info

     !>  Simple data storage at each level
     real(pfdp), allocatable :: &
          send(:),    &                 !! send buffer
          recv(:),    &                 !! recv buffer
          nodes(:),   &                 !! list of SDC nodes
          rmat(:,:),  &                 !! time restriction matrix
          tmat(:,:)                     !! time interpolation matrix

     integer, allocatable :: &
          nflags(:)                     !! sdc node flags

     !>  Solution variable storage
     class(pf_encap_t), allocatable :: &
          Q(:),     &           !! solution at sdc nodes
          pQ(:),    &           !! unknowns at sdc nodes, previous sweep
          R(:),     &           !! full residuals
          I(:),     &           !! 0 to node integrals
          Fflt(:),  &           !! functions values at sdc nodes (flat)
          tauQ(:),  &           !! fas correction in Q form
          pFflt(:), &           !! functions at sdc nodes, previous sweep (flat)
          q0,       &           !! initial condition 
          qend                  !! solution at end time

     !>  Function  storage
     class(pf_encap_t), pointer :: &
          F(:,:), &                     !! functions values at sdc nodes
          pF(:,:)                       !! functions at sdc nodes, previous sweep


     integer, allocatable :: shape(:)   !! user defined shape array
     type(pf_sdcmats_t), allocatable :: sdcmats
     logical :: allocated = .false.
  end type pf_level_t