add fast mode to many other inegrators

Author: szaghi

src/lib/foodie_integrand_object.f90

@@ -52,13 +52,13 @@ module foodie_integrand_object
     ! public methods for fast operational mode, must be overridden
     procedure, pass(self), public :: t_fast !< Time derivative, residuals, fast mode.
     procedure, pass(opr), public :: integrand_add_integrand_fast !< `+` fast operator.
-    generic, public :: add_fast => integrand_add_integrand_fast !< Overloading `add` method.
+    generic, public :: add_fast => integrand_add_integrand_fast !< Overloading `add_fast` method.
     procedure, pass(opr), public :: integrand_multiply_integrand_fast   !< `*` fast operator.
     procedure, pass(opr), public :: integrand_multiply_real_scalar_fast !< `* real_scalar` fast operator.
     generic, public :: multiply_fast => integrand_multiply_integrand_fast, &
-                                        integrand_multiply_real_scalar_fast !< Overloading `multiply` method.
+                                        integrand_multiply_real_scalar_fast !< Overloading `multiply_fast` method.
     procedure, pass(opr), public :: integrand_subtract_integrand_fast !< `-` fast operator.
-    generic, public :: subtract_fast => integrand_multiply_integrand_fast !< Overloading `subtract` method.
+    generic, public :: subtract_fast => integrand_subtract_integrand_fast !< Overloading `subtract_fast` method.
 endtype integrand_object
 
 abstract interface

src/lib/foodie_integrator_backward_differentiation_formula.f90

@@ -50,7 +50,7 @@ module foodie_integrator_backward_differentiation_formula
                                                          trim(class_name_)//'_5', &
                                                          trim(class_name_)//'_6'] !< List of supported schemes.
 
-logical, parameter :: has_fast_mode_=.false. !< Flag to check if integrator provides *fast mode* integrate.
+logical, parameter :: has_fast_mode_=.true. !< Flag to check if integrator provides *fast mode* integrate.
 
 type, extends(integrator_object) :: integrator_back_df
   !< FOODIE integrator: provide an implicit class of Backward-Differentiation-Formula multi-step schemes, from 1st to 6th order
@@ -73,6 +73,7 @@ module foodie_integrator_backward_differentiation_formula
     procedure, pass(self) :: destroy         !< Destroy the integrator.
     procedure, pass(self) :: initialize      !< Initialize (create) the integrator.
     procedure, pass(self) :: integrate       !< Integrate integrand field.
+    procedure, pass(self) :: integrate_fast  !< Integrate integrand field, fast mode.
     procedure, pass(self) :: update_previous !< Cyclic update previous time steps.
 endtype integrator_back_df
 
@@ -246,6 +247,38 @@ subroutine integrate(self, U, previous, Dt, t, iterations, autoupdate)
   if (autoupdate_) call self%update_previous(U=U, previous=previous)
   endsubroutine integrate
 
+  subroutine integrate_fast(self, U, previous, buffer, Dt, t, iterations, autoupdate)
+  !< Integrate field with BDF class scheme.
+  class(integrator_back_df),  intent(in)    :: self         !< Integrator.
+  class(integrand_object),    intent(inout) :: U            !< Field to be integrated.
+  class(integrand_object),    intent(inout) :: previous(1:) !< Previous time steps solutions of integrand field.
+  class(integrand_object),    intent(inout) :: buffer       !< Temporary buffer for doing fast operation.
+  real(R_P),                  intent(in)    :: Dt           !< Time steps.
+  real(R_P),                  intent(in)    :: t(:)         !< Times.
+  integer(I_P), optional,     intent(in)    :: iterations   !< Fixed point iterations.
+  logical,      optional,     intent(in)    :: autoupdate   !< Perform cyclic autoupdate of previous time steps.
+  integer(I_P)                              :: iterations_  !< Fixed point iterations.
+  logical                                   :: autoupdate_  !< Perform cyclic autoupdate of previous time steps, dummy var.
+  class(integrand_object), allocatable      :: delta        !< Delta RHS for fixed point iterations.
+  integer(I_P)                              :: s            !< Steps counter.
+
+  autoupdate_ = .true. ; if (present(autoupdate)) autoupdate_ = autoupdate
+  iterations_ = 1 ; if (present(iterations)) iterations_ = iterations
+  allocate(delta, mold=U)
+  call delta%multiply_fast(lhs=previous(self%steps), rhs=-self%a(self%steps))
+  do s=1, self%steps - 1
+    call buffer%multiply_fast(lhs=previous(s), rhs=-self%a(s))
+    call delta%add_fast(lhs=delta, rhs=buffer)
+  enddo
+  do s=1, iterations
+    buffer = U
+    call buffer%t_fast(t=t(self%steps) + Dt)
+    call buffer%multiply_fast(lhs=buffer, rhs=Dt * self%b)
+    call U%add_fast(lhs=delta, rhs=buffer)
+  enddo
+  if (autoupdate_) call self%update_previous(U=U, previous=previous)
+  endsubroutine integrate_fast
+
   subroutine update_previous(self, U, previous)
   !< Cyclic update previous time steps.
   class(integrator_back_df), intent(in)    :: self         !< Integrator.

src/lib/foodie_integrator_euler_explicit.f90

@@ -25,7 +25,7 @@ module foodie_integrator_euler_explicit
 character(len=99), parameter :: class_name_='euler_explicit'                !< Name of the class of schemes.
 character(len=99), parameter :: supported_schemes_(1:1)=[trim(class_name_)] !< List of supported schemes.
 
-logical, parameter :: has_fast_mode_=.false. !< Flag to check if integrator provides *fast mode* integrate.
+logical, parameter :: has_fast_mode_=.true. !< Flag to check if integrator provides *fast mode* integrate.
 
 type, extends(integrator_object) :: integrator_euler_explicit
   !< FOODIE integrator: provide explicit Euler scheme, it being 1st order accurate.
@@ -40,8 +40,9 @@ module foodie_integrator_euler_explicit
     procedure, pass(self) :: is_supported         !< Return .true. if the integrator class support the given scheme.
     procedure, pass(self) :: supported_schemes    !< Return the list of supported schemes.
     ! public methods
-    procedure, pass(self) :: destroy   !< Destroy the integrator.
-    procedure, nopass     :: integrate !< Integrate integrand field.
+    procedure, pass(self) :: destroy        !< Destroy the integrator.
+    procedure, nopass     :: integrate      !< Integrate integrand field.
+    procedure, nopass     :: integrate_fast !< Integrate integrand field, fast mode.
 endtype integrator_euler_explicit
 contains
   ! deferred methods
@@ -121,4 +122,17 @@ subroutine integrate(U, Dt, t)
 
   U = U + (U%t(t=t) * Dt)
   endsubroutine integrate
+
+  subroutine integrate_fast(U, buffer, Dt, t)
+  !< Integrate field with explicit Euler scheme, 1st order, fast mode.
+  class(integrand_object), intent(inout) :: U      !< Field to be integrated.
+  class(integrand_object), intent(inout) :: buffer !< Temporary buffer for doing fast operation.
+  real(R_P),               intent(in)    :: Dt     !< Time step.
+  real(R_P), optional,     intent(in)    :: t      !< Time.
+
+  buffer = U
+  call buffer%t_fast(t=t)
+  call buffer%multiply_fast(lhs=buffer, rhs=Dt)
+  call U%add_fast(lhs=U, rhs=buffer)
+  endsubroutine integrate_fast
 endmodule foodie_integrator_euler_explicit

src/lib/foodie_integrator_leapfrog.f90

@@ -49,7 +49,7 @@ module foodie_integrator_leapfrog
 character(len=99), parameter :: supported_schemes_(1:2)=[trim(class_name_)//'    ', &
                                                          trim(class_name_)//'_raw'] !< List of supported schemes.
 
-logical, parameter :: has_fast_mode_=.false. !< Flag to check if integrator provides *fast mode* integrate.
+logical, parameter :: has_fast_mode_=.true. !< Flag to check if integrator provides *fast mode* integrate.
 
 type, extends(integrator_object) :: integrator_leapfrog
   !< FOODIE integrator: provide an explicit class of leapfrog multi-step schemes, 2nd order accurate.
@@ -69,9 +69,10 @@ module foodie_integrator_leapfrog
     procedure, pass(self) :: is_supported         !< Return .true. if the integrator class support the given scheme.
     procedure, pass(self) :: supported_schemes    !< Return the list of supported schemes.
     ! public methods
-    procedure, pass(self) :: destroy    !< Destroy the integrator.
-    procedure, pass(self) :: initialize !< Initialize (create) the integrator.
-    procedure, pass(self) :: integrate  !< Integrate integrand field.
+    procedure, pass(self) :: destroy        !< Destroy the integrator.
+    procedure, pass(self) :: initialize     !< Initialize (create) the integrator.
+    procedure, pass(self) :: integrate      !< Integrate integrand field.
+    procedure, pass(self) :: integrate_fast !< Integrate integrand field, fast mode.
 endtype integrator_leapfrog
 
 contains
@@ -189,4 +190,34 @@ subroutine integrate(self, U, previous, Dt, t, filter)
   previous(1) = previous(2)
   previous(2) = U
   endsubroutine integrate
+
+  subroutine integrate_fast(self, U, previous, buffer, Dt, t, filter)
+  !< Integrate field with leapfrog class scheme, fast mode.
+  class(integrator_leapfrog),        intent(in)    :: self          !< Integrator.
+  class(integrand_object),           intent(inout) :: U             !< Field to be integrated.
+  class(integrand_object),           intent(inout) :: previous(1:2) !< Previous time steps solutions of integrand field.
+  class(integrand_object),           intent(inout) :: buffer        !< Temporary buffer for doing fast operation.
+  real(R_P),                         intent(in)    :: Dt            !< Time step.
+  real(R_P),                         intent(in)    :: t             !< Time.
+  class(integrand_object), optional, intent(inout) :: filter        !< Filter field displacement.
+
+  buffer = previous(2)
+  call buffer%t_fast(t=t)
+  call buffer%multiply_fast(lhs=buffer, rhs=Dt * 2._R_P)
+  call U%add_fast(lhs=previous(1), rhs=buffer)
+  if (present(filter)) then
+    call buffer%multiply_fast(lhs=previous(2), rhs=2._R_P)
+    call buffer%subtract_fast(lhs=previous(1), rhs=buffer)
+    call buffer%add_fast(lhs=buffer, rhs=U)
+    call filter%multiply_fast(lhs=buffer, rhs=self%nu * 0.5_R_P)
+
+    call buffer%multiply_fast(lhs=filter, rhs=self%alpha)
+    call previous(2)%add_fast(lhs=previous(2), rhs=buffer)
+
+    call buffer%multiply_fast(lhs=filter, rhs=self%alpha - 1._R_P)
+    call U%add_fast(lhs=U, rhs=buffer)
+  endif
+  previous(1) = previous(2)
+  previous(2) = U
+  endsubroutine integrate_fast
 endmodule foodie_integrator_leapfrog

src/lib/foodie_integrator_lmm_ssp.f90

@@ -40,7 +40,7 @@ module foodie_integrator_lmm_ssp
                                                          trim(class_name_)//'_steps_4_order_3', &
                                                          trim(class_name_)//'_steps_5_order_3'] !< List of supported schemes.
 
-logical, parameter :: has_fast_mode_=.false. !< Flag to check if integrator provides *fast mode* integrate.
+logical, parameter :: has_fast_mode_=.true. !< Flag to check if integrator provides *fast mode* integrate.
 
 type, extends(integrator_object) :: integrator_lmm_ssp
   !< FOODIE integrator: provide an explicit class of Linear Multi-step Methods (LLM) with Strong Stability Preserving property,
@@ -63,6 +63,7 @@ module foodie_integrator_lmm_ssp
     procedure, pass(self) :: destroy         !< Destroy the integrator.
     procedure, pass(self) :: initialize      !< Initialize (create) the integrator.
     procedure, pass(self) :: integrate       !< Integrate integrand field.
+    procedure, pass(self) :: integrate_fast  !< Integrate integrand field, fast mode.
     procedure, pass(self) :: update_previous !< Cyclic update previous time steps.
 endtype integrator_lmm_ssp
 
@@ -231,6 +232,35 @@ subroutine integrate(self, U, previous, Dt, t, autoupdate)
   if (autoupdate_) call self%update_previous(U=U, previous=previous)
   endsubroutine integrate
 
+  subroutine integrate_fast(self, U, previous, buffer, Dt, t, autoupdate)
+  !< Integrate field with LMM-SSP class scheme, fast mode.
+  class(integrator_lmm_ssp), intent(in)    :: self         !< Integrator.
+  class(integrand_object),   intent(inout) :: U            !< Field to be integrated.
+  class(integrand_object),   intent(inout) :: previous(1:) !< Previous time steps solutions of integrand field.
+  class(integrand_object),   intent(inout) :: buffer       !< Temporary buffer for doing fast operation.
+  real(R_P),                 intent(in)    :: Dt           !< Time steps.
+  real(R_P),                 intent(in)    :: t(:)         !< Times.
+  logical, optional,         intent(in)    :: autoupdate   !< Perform cyclic autoupdate of previous time steps.
+  logical                                  :: autoupdate_  !< Perform cyclic autoupdate of previous time steps, dummy var.
+  integer(I_P)                             :: s            !< Steps counter.
+
+  autoupdate_ = .true. ; if (present(autoupdate)) autoupdate_ = autoupdate
+  call U%multiply_fast(lhs=U, rhs=0._R_P)
+  do s=1, self%steps
+    if (self%a(s) /= 0._R_P) then
+       call buffer%multiply_fast(lhs=previous(s), rhs=self%a(s))
+       call U%add_fast(lhs=U, rhs=buffer)
+    endif
+    if (self%b(s) /= 0._R_P) then
+       buffer = previous(s)
+       call buffer%t_fast(t=t(s))
+       call buffer%multiply_fast(lhs=buffer, rhs=self%b(s))
+       call U%add_fast(lhs=U, rhs=buffer)
+    endif
+  enddo
+  if (autoupdate_) call self%update_previous(U=U, previous=previous)
+  endsubroutine integrate_fast
+
   subroutine update_previous(self, U, previous)
   !< Cyclic update previous time steps.
   class(integrator_lmm_ssp), intent(in)    :: self         !< Integrator.

src/lib/foodie_integrator_lmm_ssp_vss.f90

@@ -357,11 +357,18 @@ subroutine integrate_order_3_fast(self, U, previous, buffer, Dt, t, autoupdate)
 
   autoupdate_ = .true. ; if (present(autoupdate)) autoupdate_ = autoupdate
   omega_= omega(Dt=Dt, s=self%steps-1)
-  ! @TODO to be completed
-  U = (previous(1) * ((3._R_P * omega_ + 2._R_P) / omega_ ** 3 )) +                           &
-      (previous(self%steps) * (((omega_ + 1._R_P) ** 2) * (omega_ - 2._R_P) / omega_ ** 3)) + &
-      (previous(1)%t(t=t(1)) * (Dt(self%steps) * (omega_ + 1._R_P) / omega_ ** 2)) +          &
-      (previous(self%steps)%t(t=t(self%steps)) * (Dt(self%steps) * (omega_ + 1._R_P) ** 2 / omega_ ** 2))
+
+  call U%multiply_fast(lhs=previous(1), rhs=(3._R_P * omega_ + 2._R_P) / (omega_ ** 3))
+  call buffer%multiply_fast(lhs=previous(self%steps), rhs=(((omega_ + 1._R_P) ** 2) * (omega_ - 2._R_P) / (omega_ ** 3)))
+  call U%add_fast(lhs=U, rhs=buffer)
+  buffer = previous(1)
+  call buffer%t_fast(t=t(1))
+  call buffer%multiply_fast(lhs=buffer, rhs=Dt(self%steps) * (omega_ + 1._R_P) / (omega_ ** 2))
+  call U%add_fast(lhs=U, rhs=buffer)
+  buffer = previous(self%steps)
+  call buffer%t_fast(t=t(self%steps))
+  call buffer%multiply_fast(lhs=buffer, rhs=(Dt(self%steps) * (omega_ + 1._R_P) ** 2 / (omega_ ** 2)))
+  call U%add_fast(lhs=U, rhs=buffer)
   if (autoupdate_) call self%update_previous(U=U, previous=previous, Dt=Dt)
   endsubroutine integrate_order_3_fast