From 3b5b887aef93aaa9e391c60011a08d2173d6a447 Mon Sep 17 00:00:00 2001
From: Maren Kaluza <maren.kaluza@ufz.de>
Date: Wed, 2 Oct 2024 11:56:48 +0200
Subject: [PATCH] Updated to type opti_sim_t with flexible opti data, minimal
 eval interface

New features:

- The eval interface only has two variables: config, opti_sim(:)
    - config contains the parameters. Also for example the indices. The
      indices are helpful for parallelization: the eval would only use a
      subset of domains.
    - opti_sim is an array that will have the size nDomains in case of
      mHM. Each element of the array is then of type opti_sim and can
      have an array of different optimization variables 'opti_vars'.
- opti_vars replaces specific optimization variables like runoff, et,
  tws, etc. You can add an arbitrary variable with dimension 1, 2, 3, 4
  or 5 with an arbitrary name like "runoff", "et", "tws", etc.
    - optivars contains:
      - add: to add a variable. -- This might not work correctly if
	a second variable is added and so on.
      - has: by name, if a variable is added.
      - allocate: allocates the space for the variable. It cannot be
	done through a pointer.
      - set_pointer: sets a pointer to the allocated space of a
	variable.

Changes of plan:

- The opti_sim input variable in the eval interface was not supposed to
  be an array. But it was more convenient that way and would have added
  another complicated layer to the code otherwise.
- There was not supposed to be a procedure 'allocate' in opti_sim. But
  int turned out, associated arrays cannot be allocated through
  pointers.
---
 src/mo_cost.f90               |  14 +-
 src/mo_likelihood.f90         |  85 ++++++++---
 src/mo_opt_functions.f90      | 147 ++++++++++++-------
 src/mo_optimization_types.f90 | 262 ++++++++++++++++++++++++++++++++--
 src/mo_optimization_utils.f90 |  20 ++-
 5 files changed, 432 insertions(+), 96 deletions(-)

diff --git a/src/mo_cost.f90 b/src/mo_cost.f90
index 1f77921..d9bbf65 100644
--- a/src/mo_cost.f90
+++ b/src/mo_cost.f90
@@ -223,9 +223,9 @@ CONTAINS
   !> \brief dummy cost objective function
   FUNCTION cost_objective(parameterset, eval, arg1, arg2, arg3)
 
-    use mo_kind, only: dp
+    use mo_kind, only: dp, i4
     use mo_optimization_utils, only: eval_interface
-    use mo_optimization_types, only : variables_optidata_sim
+    use mo_optimization_types, only : opti_sim_t, config_t
 
     implicit none
 
@@ -236,11 +236,15 @@ CONTAINS
     real(dp), optional, intent(out) :: arg3
     real(dp) :: cost_objective
 
-    type(variables_optidata_sim) :: et_opti
+    type(opti_sim_t), dimension(:), pointer :: opti_sim
+    type(config_t) :: config
     REAL(DP), DIMENSION(6,2)            :: meas
     REAL(DP), DIMENSION(6)              :: calc
 
-    call eval(parameterset, et_opti)
+    config%parameters = parameterset
+    allocate(opti_sim(1))
+    call opti_sim(1)%add(name='et', dim=2_i4)
+    call eval(config, opti_sim)
 
     ! function: f(x) = ax^3 + bx^2 + cx + d
     ! measurements: (0.5,5.725), (1.0, 21.7), (1.5, 49.175), (2.0, 88.9), (2.5, 141.625), (3.0, 208.1)
@@ -254,6 +258,8 @@ CONTAINS
     ! MAE  Mean Absolute Error
     cost_objective = sum(abs( meas(:,2)-calc(:) ))/size(meas,1)
 
+    deallocate(opti_sim)
+
     RETURN
   END FUNCTION cost_objective
 
diff --git a/src/mo_likelihood.f90 b/src/mo_likelihood.f90
index 4a7b70e..74d33fd 100644
--- a/src/mo_likelihood.f90
+++ b/src/mo_likelihood.f90
@@ -47,7 +47,7 @@ CONTAINS
   ! -------------------------------
   !> \brief A Likelihood function: "real" likelihood  (sigma is an error model or given)
   function likelihood_dp(paraset, eval, stddev_in, stddev_new, likeli_new)
-    use mo_optimization_types, only : variables_optidata_sim
+    use mo_optimization_types, only : opti_sim_t, config_t
     REAL(DP), DIMENSION(:), INTENT(IN)            :: paraset          ! parameter set
     procedure(eval_interface), INTENT(IN), pointer :: eval
     REAL(DP),               INTENT(IN), optional  :: stddev_in        ! standard deviation of data
@@ -58,18 +58,26 @@ CONTAINS
 
     ! local
     REAL(DP), DIMENSION(size(meas,1))   :: errors
-    type(variables_optidata_sim) :: runoff
+    real(dp), pointer :: runoff(:, :)
+    type(opti_sim_t), dimension(:), pointer :: opti_sim
+    type(config_t) :: config
 
-    call eval(paraset, runoff)
-    errors(:) = runoff%runoff(:,1)-data()
+    config%parameters = paraset
+
+    allocate(opti_sim(1))
+    call opti_sim(1)%add(name="runoff", dim=2_i4)
+    call eval(config, opti_sim=opti_sim)
+    call opti_sim(1)%set_pointer(ptr=runoff, name="runoff")
+    errors(:) = runoff(:,1)-data()
     likelihood_dp = exp(-0.5_dp * sum( errors(:) * errors(:) / stddev_global**2 ))
+    deallocate(opti_sim)
 
   end function likelihood_dp
 
   ! -------------------------------
   !> \brief A Log-Likelihood function: "real" likelihood  (sigma is an error model or given)
   function loglikelihood_dp(paraset, eval, stddev_in, stddev_new, likeli_new)
-    use mo_optimization_types, only : variables_optidata_sim
+    use mo_optimization_types, only : opti_sim_t, config_t
     REAL(DP), DIMENSION(:), INTENT(IN)            :: paraset          ! parameter set
     procedure(eval_interface), INTENT(IN), pointer :: eval
     REAL(DP),               INTENT(IN), optional  :: stddev_in        ! standard deviation of data
@@ -80,18 +88,26 @@ CONTAINS
 
     ! local
     REAL(DP), DIMENSION(size(meas,1))   :: errors
-    type(variables_optidata_sim) :: runoff
+    real(dp), pointer :: runoff(:, :)
+    type(opti_sim_t), dimension(:), pointer :: opti_sim
+    type(config_t) :: config
+
+    config%parameters = paraset
 
-    call eval(paraset, runoff)
-    errors(:) = runoff%runoff(:,1)-data()
+    allocate(opti_sim(1))
+    call opti_sim(1)%add(name="runoff", dim=2_i4)
+    call eval(config, opti_sim=opti_sim)
+    call opti_sim(1)%set_pointer(ptr=runoff, name="runoff")
+    errors(:) = runoff(:,1)-data()
     loglikelihood_dp = -0.5_dp * sum( errors(:) * errors(:) / stddev_global**2 )
+    deallocate(opti_sim)
 
   end function loglikelihood_dp
 
   ! -------------------------------
   !> \brief A Likelihood function: "faked" likelihood (sigma is computed by obs vs model)
   function likelihood_stddev_dp(paraset, eval, stddev_in, stddev_new, likeli_new)
-    use mo_optimization_types, only : variables_optidata_sim
+    use mo_optimization_types, only : opti_sim_t, config_t
     REAL(DP), DIMENSION(:), INTENT(IN)            :: paraset          ! parameter set
     procedure(eval_interface), INTENT(IN), pointer :: eval
     REAL(DP),               INTENT(IN), optional  :: stddev_in        ! standard deviation of data
@@ -103,11 +119,19 @@ CONTAINS
     ! local
     REAL(DP), DIMENSION(size(meas,1))   :: errors
     REAL(DP)                            :: stddev_err
-    type(variables_optidata_sim) :: runoff
+    real(dp), pointer :: runoff(:, :)
+    type(opti_sim_t), dimension(:), pointer :: opti_sim
+    type(config_t) :: config
 
-    call eval(paraset, runoff)
-    errors(:) = runoff%runoff(:,1)-data()
+    config%parameters = paraset
+
+    allocate(opti_sim(1))
+    call opti_sim(1)%add(name="runoff", dim=2_i4)
+    call eval(config, opti_sim=opti_sim)
+    call opti_sim(1)%set_pointer(ptr=runoff, name="runoff")
+    errors(:) = runoff(:,1)-data()
     likelihood_stddev_dp = exp(-0.5_dp * sum( errors(:) * errors(:) / stddev_in**2 ))
+    deallocate(opti_sim)
 
     ! optional out
     stddev_err = stddev(errors)
@@ -123,7 +147,7 @@ CONTAINS
   ! -------------------------------
   !> \brief A Log-Likelihood_stddev function: "faked" likelihood (sigma is computed by obs vs model)
   function loglikelihood_stddev_dp(paraset, eval, stddev_in, stddev_new, likeli_new)
-    use mo_optimization_types, only : variables_optidata_sim
+    use mo_optimization_types, only : opti_sim_t, config_t
     REAL(DP), DIMENSION(:), INTENT(IN)            :: paraset          ! parameter set
     procedure(eval_interface), INTENT(IN), pointer :: eval
     REAL(DP),               INTENT(IN), optional  :: stddev_in        ! standard deviation of data
@@ -135,11 +159,19 @@ CONTAINS
     ! local
     REAL(DP), DIMENSION(size(meas,1))   :: errors
     REAL(DP)                            :: stddev_err
-    type(variables_optidata_sim) :: runoff
+    real(dp), pointer :: runoff(:, :)
+    type(opti_sim_t), dimension(:), pointer :: opti_sim
+    type(config_t) :: config
+
+    config%parameters = paraset
 
-    call eval(paraset, runoff)
-    errors(:) = runoff%runoff(:,1)-data()
+    allocate(opti_sim(1))
+    call opti_sim(1)%add(name="runoff", dim=2_i4)
+    call eval(config, opti_sim=opti_sim)
+    call opti_sim(1)%set_pointer(ptr=runoff, name="runoff")
+    errors(:) = runoff(:,1)-data()
     loglikelihood_stddev_dp = -0.5_dp * sum( errors(:) * errors(:) / stddev_in**2 )
+    deallocate(opti_sim)
 
     ! optional out
     stddev_err = stddev(errors)
@@ -154,21 +186,28 @@ CONTAINS
 
   ! -------------------------------
   !> \brief A Model: p1*x^2 + p2*x + p3
-  subroutine model_dp(parameterset, varsOptidataSim)
+  subroutine model_dp(config, opti_sim)
 
     use mo_kind, only: dp
-    use mo_optimization_types, only : variables_optidata_sim
+    use mo_optimization_types, only : opti_sim_t, config_t
     !! !$ USE omp_lib,    only: OMP_GET_THREAD_NUM
 
-    real(dp),    dimension(:), intent(in) :: parameterset
-    type(variables_optidata_sim), intent(inout) :: varsOptidataSim
+    type(config_t), intent(in) :: config
+    type(opti_sim_t), dimension(:), pointer, optional, intent(inout) :: opti_sim
+    real(dp), pointer :: runoff(:, :)
 
     integer(i4) :: i, n
     ! for OMP
     !! !$  integer(i4)                           :: n_threads, is_thread
 
     n = size(meas,1)
-    allocate(varsOptidataSim%runoff(n, 1))
+
+    do i = 1 , size(opti_sim)
+      if (opti_sim(i)%has('runoff')) then
+        call opti_sim(i)%allocate(name="runoff", dim1=n, dim2=1)
+        call opti_sim(i)%set_pointer(ptr=runoff, name="runoff")
+      end if
+    end do
 
     !! !$ is_thread = OMP_GET_THREAD_NUM()
     !! !$ write(*,*) 'OMP_thread: ', is_thread
@@ -178,11 +217,11 @@ CONTAINS
     !$OMP do
     do i=1, n
        !! !$ if (is_thread /= 0) write(*,*) '    OMP_thread-1: ', is_thread
-       varsOptidataSim%runoff(i,1) = parameterset(1) * meas(i,1) * meas(i,1) + parameterset(2) * meas(i,1) + parameterset(3)
+       runoff(i,1) = config%parameters(1) * meas(i,1) * meas(i,1) + config%parameters(2) * meas(i,1) + config%parameters(3)
     end do
     !$OMP end do
     !$OMP end parallel
-
+    
   end subroutine model_dp
 
   function data_dp()
diff --git a/src/mo_opt_functions.f90 b/src/mo_opt_functions.f90
index 4a74802..6a9fb66 100644
--- a/src/mo_opt_functions.f90
+++ b/src/mo_opt_functions.f90
@@ -5623,7 +5623,7 @@ CONTAINS
   function ackley_objective(parameterset, eval, arg1, arg2, arg3)
 
     use mo_constants, only: pi_dp
-    use mo_optimization_types, only : variables_optidata_sim
+    use mo_optimization_types, only : opti_sim_t, config_t
 
     implicit none
 
@@ -5639,21 +5639,27 @@ CONTAINS
     real(dp), parameter :: b = 0.2_dp
     real(dp), parameter :: c = 2.0_dp*pi_dp
     real(dp) :: s1, s2
-    type(variables_optidata_sim) :: et_opti
+    type(opti_sim_t), pointer, dimension(:) :: opti_sim
+    type(config_t) :: config
 
-    call eval(parameterset, et_opti)
+    allocate(opti_sim(1))
+    config%parameters = parameterset
+    call opti_sim(1)%add(name='et', dim=2_i4)
+    call eval(config, opti_sim)
 
     n = size(parameterset)
     s1 = sum(parameterset**2)
     s2 = sum(cos(c*parameterset))
     ackley_objective = -a * exp(-b*sqrt(1.0_dp/real(n,dp)*s1)) - exp(1.0_dp/real(n,dp)*s2) + a + exp(1.0_dp)
 
+    deallocate(opti_sim)
+
   end function ackley_objective
 
   function griewank_objective(parameterset, eval, arg1, arg2, arg3)
 
     use mo_kind, only: i4, dp
-    use mo_optimization_types, only : variables_optidata_sim
+    use mo_optimization_types, only : opti_sim_t, config_t
 
     implicit none
 
@@ -5667,9 +5673,13 @@ CONTAINS
     integer(i4) :: nopt
     integer(i4) :: j
     real(dp)    :: d, u1, u2
-    type(variables_optidata_sim) :: et_opti
+    type(opti_sim_t), pointer, dimension(:) :: opti_sim
+    type(config_t) :: config
 
-    call eval(parameterset, et_opti)
+    config%parameters = parameterset
+    allocate(opti_sim(1))
+    call opti_sim(1)%add(name='et', dim=2_i4)
+    call eval(config, opti_sim)
 
     nopt = size(parameterset)
     if (nopt .eq. 2) then
@@ -5683,69 +5693,110 @@ CONTAINS
        u2 = u2 * cos(parameterset(j)/sqrt(real(j,dp)))
     end do
     griewank_objective = u1 - u2 + 1.0_dp
+
+    deallocate(opti_sim)
     !
   end function griewank_objective
 
 
-  subroutine eval_dummy(parameterset, varsOptidataSim)
+  subroutine eval_dummy(config, opti_sim)
     use mo_kind, only : dp
-    use mo_optimization_types, only : variables_optidata_sim, optidata
+    use mo_optimization_types, only : opti_sim_t, config_t, optidata
 
     implicit none
 
-    real(dp),    dimension(:), intent(in) :: parameterset
-    type(variables_optidata_sim), intent(inout) :: varsOptidataSim
+    type(config_t), intent(in) :: config
+    type(opti_sim_t), dimension(:), pointer, optional, intent(inout) :: opti_sim
 
-    type(optidata) :: dummyData
-    integer(i4) :: i
+    real(dp), dimension(:, :), pointer :: dummyDataPtr_2d
+    real(dp), dimension(:), pointer :: dummyDataPtr_1d
+    integer(i4) :: iDomain, nDomains
 
-    allocate(dummyData%dataObs(1, 1))
-    dummyData%dataObs = 0.0_dp
+    nDomains = size(opti_sim)
 
-    if (allocated(varsOptidataSim%etOptiSim)) then
-      do i=1, size(varsOptidataSim%etOptiSim)
-        call varsOptidataSim%etOptiSim(i)%init(dummyData)
-      end do
-    end if
+    allocate(dummyDataPtr_2d(1, 1))
+    allocate(dummyDataPtr_1d(1))
+    dummyDataPtr_2d = 0.0_dp
+    dummyDataPtr_1d = 0.0_dp
 
-    if (allocated(varsOptidataSim%neutronsOptiSim)) then
-      do i=1, size(varsOptidataSim%neutronsOptiSim)
-        call varsOptidataSim%neutronsOptiSim(1)%init(dummyData)
-      end do
-    end if
+    do iDomain = 1 , nDomains
+      if (opti_sim(iDomain)%has('et')) then
+        call opti_sim(iDomain)%allocate(name="et", dim1=1, dim2=1)
+        call opti_sim(iDomain)%set_pointer(ptr=dummyDataPtr_2d, name="et")
+      end if
+    end do
 
-    if (allocated(varsOptidataSim%twsOptiSim)) then
-      do i=1, size(varsOptidataSim%twsOptiSim)
-        call varsOptidataSim%twsOptiSim(1)%init(dummyData)
-      end do
-    end if
+    do iDomain = 1 , nDomains
+      if (opti_sim(iDomain)%has('neutrons')) then
+        call opti_sim(iDomain)%allocate(name="neutrons", dim1=1, dim2=1)
+        call opti_sim(iDomain)%set_pointer(ptr=dummyDataPtr_2d, name="neutrons")
+      end if
+    end do
 
-    if (allocated(varsOptidataSim%smOptiSim)) then
-      do i=1, size(varsOptidataSim%smOptiSim)
-        call  varsOptidataSim%smOptiSim(1)%init(dummyData)
-      end do
-    end if
+    do iDomain = 1 , nDomains
+      if (opti_sim(iDomain)%has('tws')) then
+        call opti_sim(iDomain)%allocate(name="tws", dim1=1, dim2=1)
+        call opti_sim(iDomain)%set_pointer(ptr=dummyDataPtr_2d, name="tws")
+      end if
+    end do
 
-    allocate(varsOptidataSim%runoff(1, 1))
-     varsOptidataSim%runoff(:, :) = 0.0_dp
+    do iDomain = 1 , nDomains
+      if (opti_sim(iDomain)%has('sm')) then
+        call opti_sim(iDomain)%allocate(name="sm", dim1=1, dim2=1)
+        call opti_sim(iDomain)%set_pointer(ptr=dummyDataPtr_2d, name="sm")
+      end if
+    end do
 
-    allocate(varsOptidataSim%BFI(1))
-     varsOptidataSim%BFI(:) = 0.0_dp
+    ! ToDo: runoff and all other variables were handled differently: Created wether
+    ! optional or not. Why? Should it be changed?
+    do iDomain = 1 , nDomains
+      if (opti_sim(iDomain)%has('runoff')) then
+        call opti_sim(iDomain)%allocate(name="runoff", dim1=1, dim2=1)
+        call opti_sim(iDomain)%set_pointer(ptr=dummyDataPtr_2d, name="runoff")
+      end if
+    end do
 
-    allocate(varsOptidataSim%lake_level(1, 1))
-     varsOptidataSim%lake_level(:, :) = 0.0_dp
+    do iDomain = 1 , nDomains
+      if (opti_sim(iDomain)%has('BFI')) then
+        call opti_sim(iDomain)%allocate(name="BFI", dim1=1)
+        call opti_sim(iDomain)%set_pointer(ptr=dummyDataPtr_1d, name="BFI")
+      end if
+    end do
 
-    allocate(varsOptidataSim%lake_volume(1, 1))
-     varsOptidataSim%lake_volume(:, :) = 0.0_dp
+    do iDomain = 1 , nDomains
+      if (opti_sim(iDomain)%has('lake_level')) then
+        call opti_sim(iDomain)%allocate(name="lake_level", dim1=1, dim2=1)
+        call opti_sim(iDomain)%set_pointer(ptr=dummyDataPtr_2d, name="lake_level")
+      end if
+    end do
 
-    allocate(varsOptidataSim%lake_area(1, 1))
-     varsOptidataSim%lake_area(:, :) = 0.0_dp
+    do iDomain = 1 , nDomains
+      if (opti_sim(iDomain)%has('lake_volume')) then
+        call opti_sim(iDomain)%allocate(name="lake_volume", dim1=1, dim2=1)
+        call opti_sim(iDomain)%set_pointer(ptr=dummyDataPtr_2d, name="lake_volume")
+      end if
+    end do
 
-    allocate(varsOptidataSim%lake_spill(1, 1))
-     varsOptidataSim%lake_spill(:, :) = 0.0_dp
+    do iDomain = 1 , nDomains
+      if (opti_sim(iDomain)%has('lake_area')) then
+        call opti_sim(iDomain)%allocate(name="lake_area", dim1=1, dim2=1)
+        call opti_sim(iDomain)%set_pointer(ptr=dummyDataPtr_2d, name="lake_area")
+      end if
+    end do
+
+    do iDomain = 1 , nDomains
+      if (opti_sim(iDomain)%has('lake_spill')) then
+        call opti_sim(iDomain)%allocate(name="lake_spill", dim1=1, dim2=1)
+        call opti_sim(iDomain)%set_pointer(ptr=dummyDataPtr_2d, name="lake_spill")
+      end if
+    end do
 
-    allocate(varsOptidataSim%lake_outflow(1, 1))
-     varsOptidataSim%lake_outflow(:, :) = 0.0_dp
+    do iDomain = 1 , nDomains
+      if (opti_sim(iDomain)%has('lake_outflow')) then
+        call opti_sim(iDomain)%allocate(name="lake_outflow", dim1=1, dim2=1)
+        call opti_sim(iDomain)%set_pointer(ptr=dummyDataPtr_2d, name="lake_outflow")
+      end if
+    end do
 
   end subroutine eval_dummy
 
diff --git a/src/mo_optimization_types.f90 b/src/mo_optimization_types.f90
index b7db554..abe0caa 100644
--- a/src/mo_optimization_types.f90
+++ b/src/mo_optimization_types.f90
@@ -12,24 +12,68 @@ MODULE mo_optimization_types
 
   IMPLICIT NONE
 
-  public :: optidata, optidata_sim, variables_optidata_sim
+  public :: optidata, optidata_sim, config_t, opti_sim_t
 
   private
 
-  type variables_optidata_sim
-    integer(i4), public, dimension(:), allocatable     :: opti_domain_indices
-    real(dp), public,    dimension(:, :), allocatable  :: runoff     !< dim1=time dim2=gauge
-    type(optidata_sim), public, dimension(:), allocatable  :: smOptiSim         !< dim1=ncells, dim2=time
-    type(optidata_sim), public, dimension(:), allocatable  :: neutronsOptiSim   !< dim1=ncells, dim2=time
-    type(optidata_sim), public, dimension(:), allocatable  :: etOptiSim         !< dim1=ncells, dim2=time
-    type(optidata_sim), public, dimension(:), allocatable  :: twsOptiSim        !< dim1=ncells, dim2=time
-    real(dp), public, dimension(:, :), allocatable :: lake_level    !< dim1=time dim2=lake
-    real(dp), public, dimension(:, :), allocatable :: lake_volume   !< dim1=time dim2=lake
-    real(dp), public, dimension(:, :), allocatable :: lake_area     !< dim1=time dim2=lake
-    real(dp), public, dimension(:, :), allocatable :: lake_spill    !< dim1=time dim2=lake
-    real(dp), public, dimension(:, :), allocatable :: lake_outflow  !< dim1=time dim2=lake
-    real(dp), public, dimension(:), allocatable :: BFI           !< baseflow index, dim1=domainID
-  end type variables_optidata_sim
+ ! type variables_optidata_sim
+ !   integer(i4), public, dimension(:), allocatable     :: opti_domain_indices
+ !   real(dp), public,    dimension(:, :), allocatable  :: runoff     !< dim1=time dim2=gauge
+ !   type(optidata_sim), public, dimension(:), allocatable  :: smOptiSim         !< dim1=ncells, dim2=time
+ !   type(optidata_sim), public, dimension(:), allocatable  :: neutronsOptiSim   !< dim1=ncells, dim2=time
+ !   type(optidata_sim), public, dimension(:), allocatable  :: etOptiSim         !< dim1=ncells, dim2=time
+ !   type(optidata_sim), public, dimension(:), allocatable  :: twsOptiSim        !< dim1=ncells, dim2=time
+ !   real(dp), public, dimension(:, :), allocatable :: lake_level    !< dim1=time dim2=lake
+ !   real(dp), public, dimension(:, :), allocatable :: lake_volume   !< dim1=time dim2=lake
+ !   real(dp), public, dimension(:, :), allocatable :: lake_area     !< dim1=time dim2=lake
+ !   real(dp), public, dimension(:, :), allocatable :: lake_spill    !< dim1=time dim2=lake
+ !   real(dp), public, dimension(:, :), allocatable :: lake_outflow  !< dim1=time dim2=lake
+ !   real(dp), public, dimension(:), allocatable :: BFI           !< baseflow index, dim1=domainID
+ ! end type variables_optidata_sim
+
+  type config_t
+    real(dp),    dimension(:), allocatable :: parameters
+    integer(i4), dimension(:), allocatable :: opti_indices
+  end type config_t
+
+!  type opti_sim_t
+!    type(optidata_sim), target, dimension(:), allocatable :: data
+!    character(256)                                        :: name
+!  end type opti_sim_t
+  type opti_sim_t
+    type(opti_sim_single_t), dimension(:), allocatable :: opti_vars
+    contains
+    procedure :: has => opti_sim_t_has
+    procedure :: add => opti_sim_t_add
+    procedure :: opti_sim_t_allocate_1d
+    procedure :: opti_sim_t_allocate_2d
+    procedure :: opti_sim_t_allocate_3d
+    procedure :: opti_sim_t_allocate_4d
+    procedure :: opti_sim_t_allocate_5d
+    procedure :: opti_sim_t_set_pointer_1d
+    procedure :: opti_sim_t_set_pointer_2d
+    procedure :: opti_sim_t_set_pointer_3d
+    procedure :: opti_sim_t_set_pointer_4d
+    procedure :: opti_sim_t_set_pointer_5d
+    generic   :: set_pointer => opti_sim_t_set_pointer_1d, opti_sim_t_set_pointer_2d, &
+      opti_sim_t_set_pointer_3d, opti_sim_t_set_pointer_4d, opti_sim_t_set_pointer_5d
+    generic   :: allocate => opti_sim_t_allocate_1d, opti_sim_t_allocate_2d, &
+      opti_sim_t_allocate_3d, opti_sim_t_allocate_4d, opti_sim_t_allocate_5d
+  end type opti_sim_t
+
+  type opti_sim_single_t
+    real(dp), dimension(:),             allocatable :: data_1d
+    real(dp), dimension(:, :),          allocatable :: data_2d
+    real(dp), dimension(:, :, :),       allocatable :: data_3d
+    real(dp), dimension(:, :, :, :),    allocatable :: data_4d
+    real(dp), dimension(:, :, :, :, :), allocatable :: data_5d
+    character(256)                                  :: name
+    integer(i4)                                     :: dimen
+    integer(i4)                                     :: time_avg_selector = 1_i4 !< time averaging: -3 yearly, -2 monthly, -1 daily,
+                                                                                !< 0 total, n every n timestep
+    ! contains
+    ! procedure :: add => opti_sim_single_t_add
+  end type opti_sim_single_t
 
 
   !> \brief optional data, such as sm, neutrons, et, tws
@@ -64,6 +108,194 @@ MODULE mo_optimization_types
 
   contains
 
+  ! ToDo: When to charater(*) and character(256)?
+  pure logical function opti_sim_t_has(this, name)
+    class(opti_sim_t), intent(in) :: this
+    character(*), intent(in)    :: name
+
+    integer(i4) :: i
+
+    opti_sim_t_has = .false.
+
+    do i = 1, size(this%opti_vars)
+      if (trim(this%opti_vars(i)%name) == trim(name)) opti_sim_t_has = .true.
+    end do
+  end function opti_sim_t_has
+
+  subroutine opti_sim_t_add(this, name, dim, time_avg_selector)
+    class(opti_sim_t), intent(inout) :: this
+    character(*), intent(in)    :: name
+    integer(i4), intent(in) :: dim
+    integer(i4), optional, intent(in) :: time_avg_selector
+
+    type(opti_sim_single_t) :: add_data
+
+    ! ToDo: Why name in type 256 and in input var *?
+    add_data%name = name
+    add_data%dimen = dim
+    if (present(time_avg_selector)) add_data%time_avg_selector = time_avg_selector
+    ! ToDo: is the if case needed?
+    if (allocated(this%opti_vars)) then
+      this%opti_vars = [this%opti_vars, add_data]
+    else
+      allocate(this%opti_vars(1))
+      this%opti_vars(1)=add_data
+    end if
+
+  end subroutine opti_sim_t_add
+
+  subroutine opti_sim_t_allocate_1d(this, name, dim1)
+    class(opti_sim_t), target, intent(inout) :: this
+    character(*), intent(in)    :: name
+    integer(i4), intent(in)  :: dim1
+
+    integer(i4) :: i
+
+    do i = 1, size(this%opti_vars)
+      if (this%opti_vars(i)%name == name) then
+        allocate(this%opti_vars(i)%data_1d(dim1))
+      end if
+    end do
+  end subroutine opti_sim_t_allocate_1d
+
+  subroutine opti_sim_t_allocate_2d(this, name, dim1, dim2)
+    class(opti_sim_t), target, intent(inout) :: this
+    character(*), intent(in)    :: name
+    integer(i4), intent(in)  :: dim1
+    integer(i4), intent(in)  :: dim2
+
+    integer(i4) :: i
+
+    do i = 1, size(this%opti_vars)
+      if (this%opti_vars(i)%name == name) then
+        allocate(this%opti_vars(i)%data_2d(dim1, dim2))
+      end if
+    end do
+  end subroutine opti_sim_t_allocate_2d
+
+  subroutine opti_sim_t_allocate_3d(this, name, dim1, dim2, dim3)
+    class(opti_sim_t), target, intent(inout) :: this
+    character(*), intent(in)    :: name
+    integer(i4), intent(in)  :: dim1
+    integer(i4), intent(in)  :: dim2
+    integer(i4), intent(in)  :: dim3
+
+    integer(i4) :: i
+
+    do i = 1, size(this%opti_vars)
+      if (this%opti_vars(i)%name == name) then
+        allocate(this%opti_vars(i)%data_3d(dim1, dim2, dim3))
+      end if
+    end do
+  end subroutine opti_sim_t_allocate_3d
+
+  subroutine opti_sim_t_allocate_4d(this, name, dim1, dim2, dim3, dim4)
+    class(opti_sim_t), target, intent(inout) :: this
+    character(*), intent(in)    :: name
+    integer(i4), intent(in)  :: dim1
+    integer(i4), intent(in)  :: dim2
+    integer(i4), intent(in)  :: dim3
+    integer(i4), intent(in)  :: dim4
+
+    integer(i4) :: i
+
+    do i = 1, size(this%opti_vars)
+      if (this%opti_vars(i)%name == name) then
+        allocate(this%opti_vars(i)%data_4d(dim1, dim2, dim3, dim4))
+      end if
+    end do
+  end subroutine opti_sim_t_allocate_4d
+
+  subroutine opti_sim_t_allocate_5d(this, name, dim1, dim2, dim3, dim4, dim5)
+    class(opti_sim_t), target, intent(inout) :: this
+    character(*), intent(in)    :: name
+    integer(i4), intent(in)  :: dim1
+    integer(i4), intent(in)  :: dim2
+    integer(i4), intent(in)  :: dim3
+    integer(i4), intent(in)  :: dim4
+    integer(i4), intent(in)  :: dim5
+
+    integer(i4) :: i
+
+    do i = 1, size(this%opti_vars)
+      if (this%opti_vars(i)%name == name) then
+        allocate(this%opti_vars(i)%data_5d(dim1, dim2, dim3, dim4, dim5))
+      end if
+    end do
+  end subroutine opti_sim_t_allocate_5d
+
+  ! ToDo: generate with fypp
+  subroutine opti_sim_t_set_pointer_1d(this, ptr, name)
+    class(opti_sim_t), target, intent(in) :: this
+    real(dp), dimension(:), pointer, intent(inout) :: ptr
+    character(*), intent(in)    :: name
+
+    integer(i4) :: i
+
+    do i = 1, size(this%opti_vars)
+      if (this%opti_vars(i)%name == name) then
+        ptr => this%opti_vars(i)%data_1d
+      end if
+    end do
+  end subroutine opti_sim_t_set_pointer_1d
+
+  subroutine opti_sim_t_set_pointer_2d(this, ptr, name)
+    class(opti_sim_t), target, intent(in) :: this
+    real(dp), dimension(:,:), pointer :: ptr
+    character(*), intent(in)    :: name
+
+    integer(i4) :: i
+
+    do i = 1, size(this%opti_vars)
+      if (this%opti_vars(i)%name == name) then
+        ptr => this%opti_vars(i)%data_2d
+      end if
+    end do
+  end subroutine opti_sim_t_set_pointer_2d
+
+  subroutine opti_sim_t_set_pointer_3d(this, ptr, name)
+    class(opti_sim_t), target, intent(in) :: this
+    real(dp), dimension(:,:,:), pointer, intent(inout) :: ptr
+    character(*), intent(in)    :: name
+
+    integer(i4) :: i
+
+    do i = 1, size(this%opti_vars)
+      if (this%opti_vars(i)%name == name) then
+        ptr => this%opti_vars(i)%data_3d
+      end if
+    end do
+  end subroutine opti_sim_t_set_pointer_3d
+
+  subroutine opti_sim_t_set_pointer_4d(this, ptr, name)
+    class(opti_sim_t), target, intent(in) :: this
+    real(dp), dimension(:,:,:,:), pointer, intent(inout) :: ptr
+    character(*), intent(in)    :: name
+
+    integer(i4) :: i
+
+    do i = 1, size(this%opti_vars)
+      if (this%opti_vars(i)%name == name) then
+        ptr => this%opti_vars(i)%data_4d
+      end if
+    end do
+  end subroutine opti_sim_t_set_pointer_4d
+
+  subroutine opti_sim_t_set_pointer_5d(this, ptr, name)
+    class(opti_sim_t), target, intent(in) :: this
+    real(dp), dimension(:,:,:,:,:), pointer, intent(inout) :: ptr
+    character(*), intent(in)    :: name
+
+    integer(i4) :: i
+
+    do i = 1, size(this%opti_vars)
+      if (this%opti_vars(i)%name == name) then
+        ptr => this%opti_vars(i)%data_5d
+      end if
+    end do
+  end subroutine opti_sim_t_set_pointer_5d
+
+  ! ToDo: Pass only shape instead of optidataObs
   subroutine optidata_sim_init(this, optidataObs)
     class(optidata_sim), intent(inout) :: this
     type(optidata),      intent(in)    :: optidataObs
diff --git a/src/mo_optimization_utils.f90 b/src/mo_optimization_utils.f90
index 3abe0f7..827f0c3 100644
--- a/src/mo_optimization_utils.f90
+++ b/src/mo_optimization_utils.f90
@@ -11,13 +11,21 @@ module mo_optimization_utils
 
   implicit none
 
-  !> \brief Interface for evaluation routine.
+ ! !> \brief Interface for evaluation routine.
+ ! abstract interface
+ !   subroutine eval_interface(parameterset, varsOptidataSim)
+ !     use mo_kind, only : dp, i4
+ !     use mo_optimization_types, only : variables_optidata_sim
+ !     real(dp),    dimension(:), intent(in) :: parameterset
+ !     type(variables_optidata_sim), intent(inout) :: varsOptidataSim
+ !   end subroutine
+ ! end interface
+
   abstract interface
-    subroutine eval_interface(parameterset, varsOptidataSim)
-      use mo_kind, only : dp, i4
-      use mo_optimization_types, only : variables_optidata_sim
-      real(dp),    dimension(:), intent(in) :: parameterset
-      type(variables_optidata_sim), intent(inout) :: varsOptidataSim
+    subroutine eval_interface(config, opti_sim)
+      use mo_optimization_types, only : config_t, opti_sim_t
+      type(config_t),                                        intent(in)    :: config
+      type(opti_sim_t), dimension(:), pointer, optional, intent(inout) :: opti_sim
     end subroutine
   end interface
 
-- 
GitLab