Merge branch 'master' into develop

Author: szaghi

src/lib/foodie_integrator_runge_kutta_lssp.f90

@@ -248,7 +248,7 @@ subroutine initialize(self, scheme, U, stages, stop_on_fail)
       allocate(self%alpha(1:self%stages)) ; self%alpha = 0._R_P
       call self%initialize_order_s
     endselect
-    self%description_ = trim(adjustl(scheme))//'_stages_'//trim(str(self%stages))
+    self%description_ = trim(adjustl(scheme))//'_stages_'//trim(str(self%stages, no_sign=.true.))
     self%registers = self%stages
     if (present(U)) call self%allocate_integrand_members(U=U)
   else

src/tests/tester/foodie_test_integrand_ladvection.f90

@@ -55,15 +55,19 @@ module foodie_test_integrand_ladvection
    real(R_P)                               :: a=0._R_P        !< Advection coefficient.
    real(R_P), allocatable                  :: u(:)            !< Integrand (state) variable.
    class(interpolator_object), allocatable :: interpolator    !< WENO interpolator.
+   character(99)                           :: initial_state   !< Initial state.
    contains
       ! auxiliary methods
       procedure, pass(self), public :: destroy          !< Destroy field.
-      procedure, pass(self), public :: dt => compute_dt !< Compute the current time step, by means of CFL condition.
-      procedure, pass(self), public :: exact_solution   !< Return exact solution.
+      procedure, pass(self), public :: dt => compute_dt !< Compute the current time step by means of CFL condition.
+      procedure, pass(self), public ::       compute_dx !< Compute the space step by means of CFL condition.
       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 integrand.
       procedure, pass(self), public :: parse_cli      !< Initialize from command line interface.
       procedure, nopass,     public :: set_cli        !< Set command line interface.
       ! integrand_object deferred methods
@@ -117,22 +121,24 @@ pure function compute_dt(self, final_time, t) result(Dt)
    real(R_P),                   intent(in), optional :: t          !< Time.
    real(R_P)                                         :: Dt         !< Time step.
 
-   associate(a=>self%a, Ni=>self%Ni, Dx=>self%Dx, CFL=>self%CFL)
+   associate(a=>self%a, Dx=>self%Dx, CFL=>self%CFL)
       Dt = Dx * CFL / abs(a)
       if (present(t)) then
          if ((t + Dt) > final_time) Dt = final_time - t
       endif
    endassociate
    endfunction compute_dt
 
-   pure function exact_solution(self, u0, t) result(exact)
-   !< Return exact solution.
-   class(integrand_ladvection), intent(in) :: self     !< Integrand.
-   real(R_P),                   intent(in) :: u0(1:)   !< Initial state
-   real(R_P),                   intent(in) :: t        !< Time.
-   real(R_P), allocatable                  :: exact(:) !< Exact solution.
+   pure function compute_dx(self, Dt) result(Dx)
+   !< Compute the space step step by means of CFL condition.
+   class(integrand_ladvection), intent(in) :: self       !< Advection field.
+   real(R_P),                   intent(in) :: Dt         !< Time step.
+   real(R_P)                               :: Dx         !< Space step.
 
-   endfunction exact_solution
+   associate(a=>self%a, CFL=>self%CFL)
+      Dx = Dt / CFL * abs(a)
+   endassociate
+   endfunction compute_dx
 
    pure function output(self) result(state)
    !< Output the advection field state.
@@ -154,6 +160,45 @@ pure function description(self, prefix) result(desc)
    desc = prefix//'linear_advection-Ni_'//trim(strz(self%Ni, 10))
    endfunction description
 
+   pure function error(self, t, U0)
+   !< Return error.
+   class(integrand_ladvection), intent(in)           :: self     !< Integrand.
+   real(R_P),                   intent(in)           :: t        !< Time.
+   class(integrand_object),     intent(in), optional :: U0       !< Initial conditions.
+   real(R_P), allocatable                            :: error(:) !< Error.
+   integer(I_P)                                      :: i        !< Counter.
+
+   allocate(error(1:1))
+   error = 0._R_P
+   if (present(U0)) then
+      select type(U0)
+      type is(integrand_ladvection)
+         do i=1, self%Ni
+            ! error = error + (U0%u(i) - self%u(i)) ** 2
+            error = max(error, abs(U0%u(i) - self%u(i)))
+         enddo
+      endselect
+      ! error = sqrt(error)
+   endif
+   endfunction error
+
+   pure function exact_solution(self, t, U0) result(exact)
+   !< Return exact solution.
+   class(integrand_ladvection), intent(in)           :: self     !< Integrand.
+   real(R_P),                   intent(in)           :: t        !< Time.
+   class(integrand_object),     intent(in), optional :: U0       !< Initial conditions.
+   real(R_P), allocatable                            :: exact(:) !< Exact solution.
+
+   if (present(U0)) then
+      select type(U0)
+      type is(integrand_ladvection)
+         exact = U0%u(1:self%Ni)
+      endselect
+   else
+      exact = self%u(1:self%Ni) * 0._R_P
+   endif
+   endfunction exact_solution
+
    subroutine export_tecplot(self, file_name, t, scheme, close_file)
    !< Export integrand to Tecplot file.
    class(integrand_ladvection), intent(in)           :: self            !< Advection field.
@@ -186,11 +231,24 @@ subroutine export_tecplot(self, file_name, t, scheme, close_file)
    endif
    endsubroutine export_tecplot
 
+   subroutine initialize(self, Dt)
+   !< Initialize integrand.
+   class(integrand_ladvection), intent(inout) :: self !< Integrand.
+   real(R_P),                   intent(in)    :: Dt   !< Time step.
+
+   self%Dx = self%compute_dx(Dt=Dt)
+   self%Ni = nint(self%length / self%Dx)
+
+   select case(trim(adjustl(self%initial_state)))
+   case('square_wave')
+      call self%set_square_wave_initial_state
+   endselect
+   endsubroutine initialize
+
    subroutine parse_cli(self, cli)
    !< Initialize from command line interface.
-   class(integrand_ladvection),  intent(inout) :: self          !< Advection field.
-   type(command_line_interface), intent(inout) :: cli           !< Command line interface handler.
-   character(99)                               :: initial_state !< Initial state.
+   class(integrand_ladvection),  intent(inout) :: self !< Advection field.
+   type(command_line_interface), intent(inout) :: cli  !< Command line interface handler.
 
    call self%destroy
 
@@ -201,16 +259,11 @@ subroutine parse_cli(self, cli)
    call cli%get(switch='-a', val=self%a, error=cli%error) ; if (cli%error/=0) stop
    call cli%get(switch='--length', val=self%length, error=cli%error) ; if (cli%error/=0) stop
    call cli%get(switch='--Ni', val=self%Ni, error=cli%error) ; if (cli%error/=0) stop
-   call cli%get(switch='-is', val=initial_state, error=cli%error) ; if (cli%error/=0) stop
+   call cli%get(switch='-is', val=self%initial_state, error=cli%error) ; if (cli%error/=0) stop
 
    self%Ng = (self%weno_order + 1) / 2
    self%Dx = self%length / self%Ni
 
-   select case(trim(adjustl(initial_state)))
-   case('square_wave')
-      call self%set_square_wave_initial_state
-   endselect
-
    if (self%weno_order>1) call wenoof_create(interpolator_type=trim(adjustl(self%w_scheme)), &
                                              S=self%Ng,                                      &
                                              interpolator=self%interpolator,                 &
@@ -425,6 +478,7 @@ pure subroutine assign_integrand(lhs, rhs)
       lhs%CFL        = rhs%CFL
       lhs%Ni         = rhs%Ni
       lhs%Ng         = rhs%Ng
+      lhs%length     = rhs%length
       lhs%Dx         = rhs%Dx
       lhs%a          = rhs%a
       if (allocated(rhs%u)) then
@@ -439,6 +493,7 @@ pure subroutine assign_integrand(lhs, rhs)
       else
          if (allocated(lhs%interpolator)) deallocate(lhs%interpolator)
       endif
+      lhs%initial_state = rhs%initial_state
    endselect
    endsubroutine assign_integrand
 

src/tests/tester/foodie_test_integrand_oscillation.f90

@@ -50,12 +50,13 @@ module foodie_test_integrand_oscillation
    contains
       ! auxiliary methods
       procedure, pass(self), public :: amplitude_phase !< Return amplitude and phase of the oscillation.
-      procedure, pass(self), public :: exact_solution  !< Return exact solution.
-      procedure, pass(self), public :: initialize      !< Initialize integrand.
       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
@@ -99,27 +100,6 @@ function amplitude_phase(self) result(ap)
    ap(2) = atan(-self%U(1) / self%U(2))
    endfunction amplitude_phase
 
-   pure function exact_solution(self, t) result(exact)
-   !< Return exact solution.
-   class(integrand_oscillation), intent(in) :: self       !< Integrand.
-   real(R_P),                    intent(in) :: t          !< Time.
-   real(R_P)                                :: exact(1:2) !< Exact solution.
-
-   exact(1) = self%U0(1) * cos(self%f * t) - self%U0(2) * sin(self%f * t)
-   exact(2) = self%U0(1) * sin(self%f * t) + self%U0(2) * cos(self%f * t)
-   endfunction exact_solution
-
-   pure subroutine initialize(self, U0, frequency)
-   !< Initialize integrand.
-   class(integrand_oscillation), intent(inout) :: self      !< Integrand.
-   real(R_P),                    intent(in)    :: U0(1:2)   !< Initial state.
-   real(R_P),                    intent(in)    :: frequency !< Frequency of oscillation.
-
-   self%U = U0
-   self%f = frequency
-   self%U0 = U0
-   endsubroutine initialize
-
    pure function output(self) result(state)
    !< Extract integrand state field.
    class(integrand_oscillation), intent(in) :: self       !< Integrand.
@@ -140,6 +120,28 @@ pure function description(self, prefix) result(desc)
    desc = prefix//'oscillation'
    endfunction description
 
+   pure function error(self, t, U0)
+   !< Return error.
+   class(integrand_oscillation), intent(in)           :: self     !< Integrand.
+   real(R_P),                    intent(in)           :: t        !< Time.
+   class(integrand_object),      intent(in), optional :: U0       !< Initial conditions.
+   real(R_P), allocatable                             :: error(:) !< Error.
+
+   allocate(error(1:2))
+   error = abs(self%U - self%exact_solution(t=t))
+   endfunction error
+
+   pure function exact_solution(self, t, U0) result(exact)
+   !< Return exact solution.
+   class(integrand_oscillation), intent(in)           :: self     !< Integrand.
+   real(R_P),                    intent(in)           :: t        !< Time.
+   class(integrand_object),      intent(in), optional :: U0       !< Initial conditions.
+   real(R_P), allocatable                             :: exact(:) !< Exact solution.
+
+   exact = [self%U0(1) * cos(self%f * t) - self%U0(2) * sin(self%f * t), &
+            self%U0(1) * sin(self%f * t) + self%U0(2) * cos(self%f * t)]
+   endfunction exact_solution
+
    subroutine export_tecplot(self, file_name, t, scheme, close_file)
    !< Export integrand to Tecplot file.
    class(integrand_oscillation), intent(in)           :: self            !< Advection field.
@@ -172,6 +174,14 @@ subroutine export_tecplot(self, file_name, t, scheme, close_file)
    endif
    endsubroutine export_tecplot
 
+   pure subroutine initialize(self, Dt)
+   !< Initialize integrand.
+   !<
+   !< Intentionally empty, all is done in `parse_cli` method.
+   class(integrand_oscillation), intent(inout) :: self !< Integrand.
+   real(R_P),                    intent(in)    :: Dt   !< Time step.
+   endsubroutine initialize
+
    subroutine parse_cli(self, cli)
    !< Initialize from command line interface.
    class(integrand_oscillation), intent(inout) :: self !< Advection field.

src/tests/tester/foodie_test_integrand_tester_object.f90

@@ -17,7 +17,10 @@ module foodie_test_integrand_tester_object
    !< This abstract provided some auxiliary methods useful for the tester machinery.
    contains
       procedure(description_interface),    pass(self), deferred :: description    !< Return an informative description of the test.
+      procedure(error_interface),          pass(self), deferred :: error          !< Return error.
+      procedure(exact_solution_interface), pass(self), deferred :: exact_solution !< Return exact solution.
       procedure(export_tecplot_interface), pass(self), deferred :: export_tecplot !< Export integrand to Tecplot file.
+      procedure(initialize_interface),     pass(self), deferred :: initialize     !< Initialize integrand.
       procedure(parse_cli_interface),      pass(self), deferred :: parse_cli      !< Initialize from command line interface.
       procedure(set_cli_interface),        nopass,     deferred :: set_cli        !< Set command line interface.
 endtype integrand_tester_object
@@ -32,20 +35,45 @@ pure function description_interface(self, prefix) result(desc)
    character(len=:), allocatable                        :: desc   !< Description.
    endfunction description_interface
 
+   pure function error_interface(self, t, U0) result(error)
+   !< Return error.
+   import :: integrand_object, integrand_tester_object, R_P
+   class(integrand_tester_object), intent(in)           :: self     !< Integrand.
+   real(R_P),                      intent(in)           :: t        !< Time.
+   class(integrand_object),        intent(in), optional :: U0       !< Initial conditions.
+   real(R_P), allocatable                               :: error(:) !< Error.
+   endfunction error_interface
+
+   pure function exact_solution_interface(self, t, U0) result(exact)
+   !< Return exact solution.
+   import :: integrand_object, integrand_tester_object, R_P
+   class(integrand_tester_object), intent(in)           :: self     !< Integrand.
+   real(R_P),                      intent(in)           :: t        !< Time.
+   class(integrand_object),        intent(in), optional :: U0       !< Initial conditions.
+   real(R_P), allocatable                               :: exact(:) !< Exact solution.
+   endfunction exact_solution_interface
+
    subroutine export_tecplot_interface(self, file_name, t, scheme, close_file)
    !< Export integrand to Tecplot file.
    import :: integrand_tester_object, R_P
-   class(integrand_tester_object), intent(in)           :: self       !< Advection field.
+   class(integrand_tester_object), intent(in)           :: self       !< Integrand.
    character(*),                   intent(in), optional :: file_name  !< File name.
    real(R_P),                      intent(in), optional :: t          !< Time.
    character(*),                   intent(in), optional :: scheme     !< Scheme used to integrate integrand.
    logical,                        intent(in), optional :: close_file !< Flag for closing file.
    endsubroutine export_tecplot_interface
 
+   subroutine initialize_interface(self, Dt)
+   !< Initialize integrand.
+   import :: integrand_tester_object, R_P
+   class(integrand_tester_object), intent(inout) :: self !< Integrand.
+   real(R_P),                      intent(in)    :: Dt   !< Time step.
+   endsubroutine initialize_interface
+
    subroutine parse_cli_interface(self, cli)
    !< Initialize from command line interface.
    import :: command_line_interface, integrand_tester_object
-   class(integrand_tester_object), intent(inout) :: self !< Advection field.
+   class(integrand_tester_object), intent(inout) :: self !< Integrand.
    type(command_line_interface),   intent(inout) :: cli  !< Command line interface handler.
    endsubroutine parse_cli_interface