time_interpolation.F90

Go to the documentation of this file.

!! Copyright (C) 2002-2006 M. Marques, A. Castro, A. Rubio, G. Bertsch
!! Copyright (C) 2023 S. Ohlmann
!!
!! This program is free software; you can redistribute it and/or modify
!! it under the terms of the GNU General Public License as published by
!! the Free Software Foundation; either version 2, or (at your option)
!! any later version.
!!
!! This program is distributed in the hope that it will be useful,
!! but WITHOUT ANY WARRANTY; without even the implied warranty of
!! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
!! GNU General Public License for more details.
!!
!! You should have received a copy of the GNU General Public License
!! along with this program; if not, write to the Free Software
!! Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
!! 02110-1301, USA.
!!
 
#include "global.h"
 
module time_interpolation_oct_m
  use debug_oct_m
  use global_oct_m
  use lalg_basic_oct_m
  use math_oct_m
  use mesh_oct_m
  use messages_oct_m
  use namespace_oct_m
  use profiling_oct_m
  use restart_oct_m
  use space_oct_m
 
  implicit none
 
  private
  public :: &
    time_interpolation_t
 
  type time_interpolation_t
    private
    complex(real64), allocatable :: zfield(:, :, :)
    real(real64), allocatable :: dfield(:, :, :)
    real(real64), allocatable :: times(:)
    integer :: max_depth
    integer :: depth
    integer :: np
    integer :: dim
    logical :: cmplx
    character(len=MAX_PATH_LEN) :: label
  contains
    procedure :: dtime_interpolation_add_time, ztime_interpolation_add_time
    generic :: add_time => dtime_interpolation_add_time, ztime_interpolation_add_time
    procedure :: dtime_interpolation_interpolate, ztime_interpolation_interpolate
    generic :: interpolate => dtime_interpolation_interpolate, ztime_interpolation_interpolate
    procedure :: read_restart => time_interpolation_read_restart
    procedure :: write_restart => time_interpolation_write_restart
    final :: time_interpolation_finalize
  end type time_interpolation_t
 
  interface time_interpolation_t
    procedure time_interpolation_constructor
  end interface time_interpolation_t
 
contains
  function time_interpolation_constructor(np, dim, depth, cmplx, label) result(this)
    integer,                     intent(in) :: np
    integer,                     intent(in) :: dim
    integer,                     intent(in) :: depth
    logical,                     intent(in) :: cmplx
    character(len=*),            intent(in) :: label
    class(time_interpolation_t), pointer :: this
 
    push_sub(time_interpolation_constructor)
 
    safe_allocate(this)
 
    this%np = np
    this%dim = dim
    this%max_depth = depth
    this%depth = 0
    this%cmplx = cmplx
    this%label = trim(label)
 
    if (this%cmplx) then
      safe_allocate(this%zfield(1:np, 1:dim, 1:depth))
    else
      safe_allocate(this%dfield(1:np, 1:dim, 1:depth))
    end if
    safe_allocate(this%times(1:depth))
 
    pop_sub(time_interpolation_constructor)
  end function time_interpolation_constructor
 
  subroutine time_interpolation_finalize(this)
    type(time_interpolation_t), intent(inout) :: this
 
    push_sub(time_interpolation_finalize)
 
    if (this%cmplx) then
      safe_deallocate_a(this%zfield)
    else
      safe_deallocate_a(this%dfield)
    end if
    safe_deallocate_a(this%times)
 
    pop_sub(time_interpolation_finalize)
  end subroutine time_interpolation_finalize
 
  subroutine time_interpolation_write_restart(this, mesh, space, restart, err)
    class(time_interpolation_t), intent(in)    :: this
    class(mesh_t),                intent(in)    :: mesh
    class(space_t),               intent(in)    :: space
    type(restart_t),              intent(in)    :: restart
    integer,                      intent(out)   :: err
 
    integer :: itime, idim, err_restart, iunit
    character(len=MAX_PATH_LEN) :: filename, lines(1)
 
    push_sub(time_interpolation_write_restart)
 
    ! write fields
    err = 0
    do itime = 1, this%depth
      do idim = 1, this%dim
        write(filename, '(a1,i2.2,a1,i3.3)') '_', itime, '_',  idim
        filename = "field_" // trim(this%label) // trim(filename)
        if (this%cmplx) then
          call restart%write_mesh_function(space, filename, mesh, &
            this%zfield(1:this%np, idim, itime), err_restart)
        else
          call restart%write_mesh_function(space, filename, mesh, &
            this%dfield(1:this%np, idim, itime), err_restart)
        end if
        if (err_restart /= 0) err = err + 1
      end do
    end do
 
    ! write times
    call restart%write_binary("field_times_"//trim(this%label), this%depth, this%times(1:this%depth), err_restart)
    if (err_restart /= 0) err = err + 1
 
    ! write depth
    iunit = restart%open("field_"//trim(this%label))
    write(lines(1), '(i2.2)') this%depth
    call restart%write(iunit, lines, 1, err_restart)
    if (err_restart /= 0) err = err + 1
    call restart%close(iunit)
 
    pop_sub(time_interpolation_write_restart)
  end subroutine time_interpolation_write_restart
 
  subroutine time_interpolation_read_restart(this, mesh, space, restart, err)
    class(time_interpolation_t),  intent(inout) :: this
    class(mesh_t),                intent(in)    :: mesh
    class(space_t),               intent(in)    :: space
    class(restart_t),             intent(in)    :: restart
    integer,                      intent(out)   :: err
 
    integer :: itime, idim, err_restart, iunit
    character(len=MAX_PATH_LEN) :: filename, lines(1)
 
    push_sub(time_interpolation_read_restart)
 
    ! read depth
    iunit = restart%open("field_"//trim(this%label))
    call restart%read(iunit, lines, 1, err_restart)
    if (err_restart /= 0) then
      err = err + 1
      pop_sub(time_interpolation_read_restart)
      return
    else
      read(lines(1), '(i2.2)') this%depth
    end if
    call restart%close(iunit)
 
    ! read fields
    err = 0
    do itime = 1, this%depth
      do idim = 1, this%dim
        write(filename, '(a1,i2.2,a1,i3.3)') '_', itime, '_',  idim
        filename = "field_" // trim(this%label) // trim(filename)
        if (this%cmplx) then
          call restart%read_mesh_function(space, filename, mesh, &
            this%zfield(1:this%np, idim, itime), err_restart)
        else
          call restart%read_mesh_function(space, filename, mesh, &
            this%dfield(1:this%np, idim, itime), err_restart)
        end if
        if (err_restart /= 0) err = err + 1
      end do
    end do
 
    ! read times
    call restart%read_binary("field_times_"//trim(this%label), this%depth, this%times(1:this%depth), err_restart)
    if (err_restart /= 0) err = err + 1
 
    pop_sub(time_interpolation_read_restart)
  end subroutine time_interpolation_read_restart
 
#include "real.F90"
#include "time_interpolation_inc.F90"
#include "undef.F90"
 
#include "complex.F90"
#include "time_interpolation_inc.F90"
#include "undef.F90"
 
end module time_interpolation_oct_m
 
!! Local Variables:
!! mode: f90
!! coding: utf-8
!! End: