source: palm/trunk/UTIL/chemistry/gasphase_preproc/kpp4palm/bin/kpp4palm.sh @ 4797

#!/bin/bash

# kpp4palm - script for creating gasphase module

#------------------------------------------------------------------------------#
# This file is part of the PALM model system.
#
# PALM 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 3 of the License, or (at your option) any later version.
#
# PALM 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
# PALM. If not, see <http://www.gnu.org/licenses/>.
#
# Copyright 2017-2020  Klaus Ketelsen and MPI-CH (April 2007)
# Copyright 2017-2020  Karlsruhe Institute of Technology
# Copyright 2017-2020  Leibniz Universitaet Hannover
#------------------------------------------------------------------------------#
# Nov. 2016: Initial Version of KPP chemistry convertor adapted for PALM 
# by Klaus Ketelsen
#
# This code is a modified version of KP4 (Jöckel, P., Kerkweg, A., Pozzer, A., 
# Sander, R., Tost, H., Riede, H., Baumgaertner, A., Gromov, S., and Kern, B.,
# 2010: Development cycle 2 of the Modular Earth Submodel System (MESSy2), 
# Geosci. Model Dev., 3, 717-752, https://doi.org/10.5194/gmd-3-717-2010).
# KP4 is part of the Modular Earth Submodel System (MESSy), which is is
# available under the  GNU General Public License (GPL).
#
#------------------------------------------------------------------------------#
#
# Current revisions:
# ------------------
# 
# 
# Former revisions:
# -----------------
# $Id: kpp4palm.sh 4481 2020-03-31 18:55:54Z raasch $ 
# Changed default mechanism from smog to phstatp
# 
# 
# 4451 2020-03-10 07:25:32Z raasch
# rename kpp4palm.ksh to kpp4palm.sh, convert to bash
# 
# 3800 2019-03-15 16:40:25Z forkel
# Removed edit of phot(nphot) for version edit phot(nphot), now done in fortran_file.C (15.03.2019, forkel)
# Editing with sed to add dummy statements toavoud unused variables (15.03.2019, forkel)
# 
# 3780 2019-03-05 11:19:45Z forkel
# forkel   05. March 2019: cs_mech and set_cm
# forkel   30. Oktober 2018: Integrating contents of kp4_compress into chem_gasphase_mod.f90
# forkel   25. September 2018: Added cat for $MECH to pass mechanism name to kpp4palm
# ketelsen 18. September 2018: Added cat for '#INLINE F90_GLOBAL' 
# (moved here from mechanisms/def_MECH/chem_gasphase.kpp
#
# forkel: 14. September  2018: WCOPY removed
# ketelsen: July 2018: Adaptations for vektor mode
# forkel June 2018: re-established original case of subroutine names
# forkel May 2018: additional copying of chem_gasphase_mod.f90 into $DEFDIR
# forkel 20.04.2018: removed  wlamch and wlamch_add from $KPP_SUBROUTINE_LIST
#                    (epsilon(one) is used now)
# forkel March 2017
# Re-introduced relative path for KPP_HOME
# Subroutine list adapted to lowercase subroutine names
# Added arr2, removed update_sun and k_3rd from subroutine list
# Renamed output file to chem_gasphase_mod
# Renamed this file from kp4/ksh to kpp4kpp.ksh
# changed location of def_mechanism directories to gasphase_preproc/mechanisms
# 
# 
# 2718 2018-01-02 08:49:38Z maronga
# Initial revision
# 
# 
##########################################################################
#
#
# Other notes:
# ------------# 
# Re-introduced relative path for KPP_HOME
# Subroutine list adapted to lowercase subroutine names
# Added arr2, removed update_sun and k_3rd from subroutine list
# Renamed output file to chem_gasphase_mod
# Renamed this fikle from kp4/ksh to kpp4kpp.ksh
# changed location of def_mechanism directories to GASPHASE_PREPROC/mechanisms
#
# Nov. 2016: Initial Version of KPP chemistry convertor by Klaus Ketelsen
#
#

set -eu


########################### User SetUp ####################################

export KPP_HOME=`pwd`/kpp
export KPP=$KPP_HOME/bin/kpp

BASE=`pwd`/kpp4palm

########################## End User Setup ################################

WORK=tmp_kpp4palm

# Default

MECH=phstatp
OUTDIR=`pwd`/../../../SOURCE
OUTFILE=chem_gasphase_mod
DEFDIR=`pwd`/mechanisms/def_phastatp
PREFIX=chem_gasphase_mod
MODE="scalar"
VLEN=1
KEEP="NO"
UPDT="NO"
DE_INDEX=0
DE_INDEX_FAST="YES"

export KPP_SOLVER=Rosenbrock

# get Command line option

while  getopts :m:i:fkup:o:s:vl:w:  c     # get options
do case $c in
      m)   MECH=$OPTARG;;            # mechanism

      i)   DE_INDEX=$OPTARG;;        # if set, deindexing

      f)   DE_INDEX_FAST="YES";;     # if set, fast deindexing

      k)   KEEP="YES";;              # keep Working directory

      o)   OUTDIR=$OPTARG;;          # Output directory of Generated Code = 
                                     # SOURCE directory (do not change)

      p)   PREFIX=$OPTARG;;          # Name Prefix (chem_gasphase_mod, do not change)

      s)   KPP_SOLVER=$OPTARG;;      # Chosen solver (only Rosebrock solvers work for vector mode)

      u)   UPDT="YES";;              # update mechanisms/def_$MECH/chem_gasphase_mod.f90

      v)   MODE="vector";;           # Set to vector Mode

      l)   VLEN=$OPTARG;;            # Set vector length

      w)   WORK=$OPTARG;;            # Working directory

      \?)  print ${0##*/} "unknown option:" $OPTARG
           print "USAGE: ${0##*/} [ -m dir -e -k -u -o dir -p name -s solver -v -l length -w dir ] "
           exit 1;;
   esac
done

echo MECHANISM = $MECH
echo DE_INDEX = $DE_INDEX
echo KEEP = $KEEP
echo UPDT = $UPDT
echo MODE = $MODE
echo VLEN = $VLEN

DEF_PREFIX=${PREFIX}.kpp
DEFDIR=`pwd`/mechanisms/def_$MECH
echo DEFDIR = $DEFDIR

# Create or clean working directory

MY_PWD=`pwd`
mkdir -p $WORK
rm -rf $WORK/*
cd $WORK

# kpp dependend, may be changed

KPP_FILE_LIST="Initialize Integrator LinearAlgebra Jacobian Function Rates Util"
if [[ $MODE = "vector" ]]
 then
KPP_FILE_LIST="$KPP_FILE_LIST kp4_compress_subroutines"
fi

KPP_SUBROUTINE_LIST="Initialize"
KPP_SUBROUTINE_LIST="$KPP_SUBROUTINE_LIST INTEGRATE Fun"
KPP_SUBROUTINE_LIST="$KPP_SUBROUTINE_LIST KppSolve KppDecomp"
KPP_SUBROUTINE_LIST="$KPP_SUBROUTINE_LIST Jac_SP k_arr "
KPP_SUBROUTINE_LIST="$KPP_SUBROUTINE_LIST Update_RCONST ARR2"
KPP_SUBROUTINE_LIST="$KPP_SUBROUTINE_LIST initialize_kpp_ctrl error_output"

# if [[ $MODE = "vector" && $KPP_SOLVER = "ROS2" ]]
# then
#   cp $BASE/templates/${KPP_SOLVER}_vec.f90 ${KPP_SOLVER}.f90    # get vector Solver 
# else
# #  KPP_SUBROUTINE_LIST="$KPP_SUBROUTINE_LIST FunTemplate JacTemplate Update_SUN "
#   KPP_SUBROUTINE_LIST="$KPP_SUBROUTINE_LIST WCOPY WSCAL WAXPY"
#   if [[ $MODE = "vector" ]]
#   then
#     cp $BASE/templates/${KPP_SOLVER}_vec.f90 ${KPP_SOLVER}.f90  # get vector Solver 
#   else
#     KPP_SUBROUTINE_LIST="$KPP_SUBROUTINE_LIST Rosenbrock  FunTemplate JacTemplate Update_SUN"
#   fi
# fi
 if [[ $MODE = "vector" ]]
 then
   # get vector Solver 
   cp $BASE/templates/${KPP_SOLVER}_vec.f90 ${KPP_SOLVER}.f90
   cp $BASE/templates/kp4_compress_header ${PREFIX}_kp4_compress_header.f90
   cp $BASE/templates/kp4_compress_subroutines ${PREFIX}_kp4_compress_subroutines.f90
fi

# Interface ignore list
KPP_INTERFACE_IGNORE=" "

echo " "
echo KPP_SOLVER $KPP_SOLVER
echo " "

case $KPP_SOLVER in
    ROS2) ;;

    Rosenbrock)   
      KPP_SUBROUTINE_LIST="$KPP_SUBROUTINE_LIST WSCAL WAXPY"
      if [[ $MODE != "vector" ]]
      then
         KPP_SUBROUTINE_LIST="$KPP_SUBROUTINE_LIST Rosenbrock  FunTemplate JacTemplate"
        KPP_INTERFACE_IGNORE="WAXPY"

      else 
        KPP_SUBROUTINE_LIST="$KPP_SUBROUTINE_LIST FunTemplate JacTemplate"
        KPP_SUBROUTINE_LIST="$KPP_SUBROUTINE_LIST kco_initialize kco_compress kco_finalize"
      fi;;

    rosenbrock_mz)
      KPP_SUBROUTINE_LIST="$KPP_SUBROUTINE_LIST WCOPY WSCAL WAXPY"
      KPP_SUBROUTINE_LIST="$KPP_SUBROUTINE_LIST Rosenbrock  FunTemplate JacTemplate Update_SUN";;

    rosenbrock)
      KPP_SUBROUTINE_LIST="$KPP_SUBROUTINE_LIST WCOPY WSCAL WAXPY"
      KPP_SUBROUTINE_LIST="$KPP_SUBROUTINE_LIST Rosenbrock  FunTemplate JacTemplate";;

    kpp_lsode)
      KPP_SUBROUTINE_LIST="$KPP_SUBROUTINE_LIST WCOPY WSCAL WAXPY"
      KPP_SUBROUTINE_LIST="$KPP_SUBROUTINE_LIST KppLsode DLSODE JAC_CHEM FUN_CHEM"
      KPP_INTERFACE_IGNORE="$KPP_INTERFACE_IGNORE JAC_CHEM KppDecomp KppSolve";;

    kpp_radau5)
      KPP_SUBROUTINE_LIST="$KPP_SUBROUTINE_LIST WCOPY WSCAL WAXPY FUN_CHEM JAC_CHEM SET2ZERO"
      KPP_SUBROUTINE_LIST="$KPP_SUBROUTINE_LIST RADAU5 Update_SUN"
      KPP_SUBROUTINE_LIST="$KPP_SUBROUTINE_LIST KppSolveCmplx KppDecompCmplx";;

    kpp_sdirk)
       KPP_SUBROUTINE_LIST="$KPP_SUBROUTINE_LIST WCOPY WSCAL WAXPY"
       KPP_SUBROUTINE_LIST="$KPP_SUBROUTINE_LIST SDIRK JAC_CHEM SET2ZERO FUN_CHEM"
       KPP_INTERFACE_IGNORE="$KPP_INTERFACE_IGNORE Set2zero SET2ZERO FUN_CHEM";;

    kpp_seulex)
       KPP_SUBROUTINE_LIST="$KPP_SUBROUTINE_LIST WCOPY WSCAL WAXPY"
       KPP_SUBROUTINE_LIST="$KPP_SUBROUTINE_LIST ATMSEULEX"
       KPP_SUBROUTINE_LIST="$KPP_SUBROUTINE_LIST SEULEX_ErrorMsg SEULEX_Integrator FUN_CHEM JAC_CHEM SEUL"
       KPP_INTERFACE_IGNORE="$KPP_INTERFACE_IGNORE SEULEX_Integrator SDIRK FUN_CHEM SEUL";;

   \?)  print "SORRY ONLY ROSENBROCK METHODS WORK AT THE MOMENT:" $KPP_SOLVER
        exit 1;;
esac
#mz-ak-20070509+

KPP_INCLUDE_LIST="Parameters Global JacobianSP Monitor"
if [[ $MODE = "vector" ]]
 then
KPP_INCLUDE_LIST="$KPP_INCLUDE_LIST kp4_compress_header"
fi

#Get definition Files

cp $DEFDIR/*.eqn         .
cp $DEFDIR/*.spc         .
cp $DEFDIR/${PREFIX}.kpp     .

# Global variable are defined here 
# This has the advantage that it is not necessary to include these variables in all .kpp definition files

cat  >> ${PREFIX}.kpp  <<  EOF
#INLINE F90_GLOBAL
! QVAP - Water vapor
  REAL(kind=dp) :: QVAP
! FAKT - Conversion factor
  REAL(kind=dp) :: FAKT

! CS_MECH for check of mechanism name with namelist
  CHARACTER(LEN=30) :: CS_MECH
#ENDINLINE
EOF

# Store mechanism name in file mech_list
cat  >> mech_list  <<  EOF
!   Mechanism: $MECH
!
EOF

# Store mechanism name for cs_mech
cat  >> set_cm  <<  EOF

! Set cs_mech for check with mechanism name from namelist
    cs_mech = '$MECH'
EOF

# Run kpp

$KPP $DEF_PREFIX

# Get templates for C++ program

cp $BASE/templates/module_header* .           # Use fixed Module_header
cp $BASE/templates/initialize_kpp_ctrl_template.f90 .  # CTRL kpp time stepping

# file with subroutine list for c++ program create_kpp_module

for i in $KPP_FILE_LIST
do
  echo ${PREFIX}_${i} >> file_list
done
echo initialize_kpp_ctrl_template >> file_list

# file with subroutine list for c++ program create_kpp_module

for i in $KPP_SUBROUTINE_LIST
do
  echo $i >> subroutine_list
done

# file with include list for c++ program create_kpp_module

for i in $KPP_INCLUDE_LIST
do
  echo ${PREFIX}_${i} >> include_list
done

touch interface_ignore_list
for i in $KPP_INTERFACE_IGNORE
do
  echo $i >> interface_ignore_list
done

echo start kpp4palm.exe with arguments
echo $PREFIX $MODE $VLEN $DE_INDEX $DE_INDEX_FAST

$BASE/bin/kpp4palm.exe $PREFIX $MODE $VLEN $DE_INDEX $DE_INDEX_FAST

# Add dummy statements in order to prevent warnings due to unused variables
#
sed -i -e '/cfactor =/a !  ' kk_kpp.f90
sed -i -e '/cfactor =/a BLANKS  IF ( lu_crow(1) == 1  .OR.  lu_icol(1) == 1  .OR.  lu_irow(1) == 1 )  CONTINUE ' kk_kpp.f90
sed -i -e '/cfactor =/a ! Following line is just to avoid compiler message about unused variables' kk_kpp.f90
sed -i -e '/cfactor =/a !  ' kk_kpp.f90

if [[ $MODE = "vector" ]]
then
sed -i -e '/! Computation of equation rates/i ! Following line is just to avoid compiler message about unused variables' kk_kpp.f90
sed -i -e '/! Computation of equation rates/i BLANKS  IF ( f(vl,nfix) > 0.0_dp )  CONTINUE' kk_kpp.f90
sed -i -e '/! Computation of equation rates/i !  ' kk_kpp.f90
else
sed -i -e '/! Computation of equation rates/i ! Following line is just to avoid compiler message about unused variables' kk_kpp.f90
sed -i -e '/! Computation of equation rates/i BLANKS  IF ( f(nfix) > 0.0_dp )  CONTINUE' kk_kpp.f90
sed -i -e '/! Computation of equation rates/i !  ' kk_kpp.f90
fi

if [[ $MODE = "vector" ]]
then
sed -i -e '/REAL(kind=dp) :: b/a  BLANKS  IF ( f(vl,nfix) > 0.0_dp )  CONTINUE' kk_kpp.f90
sed -i -e '/REAL(kind=dp) :: b/a ! Following line is just to avoid compiler message about unused variables' kk_kpp.f90
sed -i -e '/REAL(kind=dp):: b/a !' kk_kpp.f90
else
sed -i -e '/REAL(kind=dp):: b/a BLANKS  IF ( f(nfix) > 0.0_dp )  CONTINUE' kk_kpp.f90
sed -i -e '/REAL(kind=dp):: b/a ! Following line is just to avoid compiler message about unused variables' kk_kpp.f90
sed -i -e '/REAL(kind=dp):: b/a !' kk_kpp.f90
fi

sed -i -e '/one=1.0_dp/a BLANKS  IF ( incx == 0 )  CONTINUE' kk_kpp.f90
sed -i -e '/one=1.0_dp/a ! Following line is just to avoid compiler message about unused variables' kk_kpp.f90
sed -i -e '/one=1.0_dp/a !  ' kk_kpp.f90

sed -i -e '/IF (alpha .eq. zero)RETURN/i !  ' kk_kpp.f90
sed -i -e '/IF (alpha .eq. zero)RETURN/i ! Following line is just to avoid compiler message about unused variables' kk_kpp.f90
sed -i -e '/IF (alpha .eq. zero)RETURN/i BLANKS  IF ( incx == 0  .OR.  incy == 0 )  CONTINUE' kk_kpp.f90

sed -i -e '/INTENT(INOUT):: b(n)/a BLANKS  IF ( pivot(1) == 0 )  CONTINUE' kk_kpp.f90
sed -i -e '/INTENT(INOUT):: b(n)/a ! Following line is just to avoid compiler message about unused variables' kk_kpp.f90
sed -i -e '/INTENT(INOUT):: b(n)/a !  ' kk_kpp.f90

sed -i -e '1,$s/BLANKS /  /  ' kk_kpp.f90

if [[ -e $OUTDIR/${OUTFILE}.f90 ]] 
then 
 mv $OUTDIR/${OUTFILE}.f90 $OUTDIR/${OUTFILE}.f90.sav
fi
cp -p kk_kpp.f90    $OUTDIR/${OUTFILE}.f90
echo " "
echo "Write kpp module -- > " $OUTDIR/${OUTFILE}.f90

if [[ $UPDT = "YES" ]]
then
cp -p kk_kpp.f90    $DEFDIR/${OUTFILE}.f90
echo " "
echo "Write kpp module -- > " $DEFDIR/${OUTFILE}.f90
fi

if [[ $KEEP = "NO" ]]
then
  cd  $MY_PWD
  rm -rf $WORK
fi
exit