v1.2.1 release

CHANGELOG.md

@@ -1,5 +1,36 @@
 # CHANGELOG
 
+## v1.2.1, 2025-01-15
+
+- Added an option to remove the logLei discontinuity in NRL formulary fit
+- Added more output to BDE integrator
+- Added more unary transforms
+- Added key for controlling BDE integrator consolidation
+- Added 0-length timestep at the start of all simulations to fix extractor manipulator bug
+- Added option to ignore Jacobian in diffusion stencil
+
+### Breaking Changes
+
+- N/A
+
+### New Features
+
+- New standard textbook option to remoge logLei discontinuity in NRL formularly fit by moving branch cutoff to e^2 * Z^2 eV
+- BDE integrator now outputs the current number of substeps and progress through those
+- New "step" unary transform. 1 if the node evaluates > 0, 0 otherwise.
+- New "filterWithin" unary transform. 0 if the value of the node is not within the first two passed realParams, otherwise no effect.
+- New "BDEConsolidationInterval" key to set the number of timesteps for the integrator to attempt consolidation after. Defaults to 50 (old hardcoded value).
+- Some generators now explicitly output information about terms as they are generated
+- "ignoreJacobian" now valid option for diffusion stencil JSON interface
+
+### Bug Fixes
+
+- Fixed extractor manipulator bug with CRM models by performing a 0 length timestep at the start of all simulations before any output
+- Added a call to update all modelbound data before the 0-length timestep to avoid some edge-case model failure
+- Fixed bug where unsorted ingoing states in some transitions lead to the wrong required density variables being identified in some generators
+- Fixed bug with integrators that allow time evolution incrementing time twice when there is only one integration step
+- Fixed bug where composite integrator global timestep's requested size is overwritten if a timestep controller is used, leading to performance issues if output-driven timestepping is used.
+
 ## v1.2.0, 2024-09-19
 
 - Added CVODE integrator as an option

docker/Dockerfile_mac.txt

@@ -0,0 +1,142 @@
+FROM --platform=linux/arm64 ubuntu:latest
+LABEL maintainer="vera.oberhauser@ukaea.uk"
+LABEL version="1.2.1"
+LABEL description="Docker container for ReMKiT1D CI"
+
+# Disable Prompt During Packages Installation
+ARG DEBIAN_FRONTEND=noninteractive
+
+RUN apt update && apt install software-properties-common -y
+
+RUN add-apt-repository -y ppa:ubuntu-toolchain-r/test
+
+# Update and Installs Required Packages for ReMKiT1D and other Packages
+RUN apt update \
+        && apt install -y \
+        gfortran-11 \
+        g++-11 \
+        python3 \
+        python3-dev \
+        pip \
+        libblas-dev \
+        libhdf5-dev \
+        liblapack-dev \
+        m4 \
+        wget \
+        git \
+        gcovr \ 
+        gcc-arm-linux-gnueabihf \
+        build-essential \
+        pkg-config
+
+# Packages need to be Compiled and Installed Using Fortran-11
+# CMake Commands State the Required Compilers to use
+ENV FC=gfortran-11
+ENV CC=gcc-11
+ENV CXX=g++-11
+
+WORKDIR /home
+
+# ReMKiT1D Requires CMake-3.18.0 and pFUnit-4.2.2 has some Incompatibilities with Newer Versions of CMake (v3.25.0-rc1)
+RUN wget https://github.com/Kitware/CMake/releases/download/v3.18.0/cmake-3.18.0.tar.gz
+RUN tar -zvxf cmake-3.18.0.tar.gz
+
+WORKDIR /home/cmake-3.18.0
+
+RUN ./bootstrap -- -DCMAKE_USE_OPENSSL=OFF
+
+RUN make
+RUN make install
+WORKDIR /home
+
+# ReMKiT1D Requires MPI and mpich-3.4.2 has been Tested and is used in the main.yml Workflow Script
+RUN wget -O mpich.tar.gz https://www.mpich.org/static/downloads/3.4.2/mpich-3.4.2.tar.gz
+RUN tar xfz mpich.tar.gz
+
+RUN mkdir mpich-install
+WORKDIR /home/mpich-3.4.2
+
+# ReMKiT1D Requires gfortran-11, gcc-11 and g++-11. CMake Installs the Package to the mpich-install Folder
+RUN ./configure CC=gcc-11 CXX=g++-11 FC=gfortran-11 --prefix=/home/mpich-install  FFLAGS=-fallow-argument-mismatch --with-device=ch3
+
+RUN make
+RUN make install
+WORKDIR /home
+
+# Include the mpich-install PATH for Compilation of ReMKiT1D
+ENV PATH=/home/mpich-install/bin:$PATH
+
+# PETSc-3.17 is the Newest Release that has been Tested on ReMKiT1D (3.18 not Compatible)
+RUN git clone -b release-3.17 https://gitlab.com/petsc/petsc.git petsc
+
+WORKDIR petsc/
+
+#PETSc Requires MPI and PETSc is Told Where MPI (mpich-4.2.2) is Installed
+RUN ./configure --with-mpi-dir=/home/mpich-install --download-hypre=1 --download-fblaslapack=1 --with-debugging=0 COPTFLAGS=-O3 CXXOPTFLAGS=-O3 FOPTFLAGS=-O3
+RUN make PETSC_DIR=/home/petsc PETSC_ARCH=arch-linux-c-opt all check
+
+WORKDIR /home
+
+# pFUnit-4.2.2 is Installed for the use of Unit Testing
+RUN git clone -b v4.2.2 https://github.com/Goddard-Fortran-Ecosystem/pFUnit.git
+WORKDIR pFUnit
+RUN mkdir build
+WORKDIR build
+
+# ReMKiT1D Requires gfortran-11, gcc-11 and g++-11. CMake Installs the Package to the 'installs' Folder
+RUN cmake .. -DSKIP_OPENMP=yes -DCMAKE_Fortran_COMPILER=gfortran-11 -DCMAKE_C_COMPILER=gcc-11 -DCMAKE_CXX_COMPILER=g++-11 -DCMAKE_INSTALL_PREFIX=/home/installs/pFUnit
+RUN make install
+WORKDIR /home
+
+# ReMKiT1D is Currently using json-fortran-8.2.5
+RUN git clone -b 8.2.5 https://github.com/jacobwilliams/json-fortran.git
+
+WORKDIR json-fortran
+RUN mkdir build
+WORKDIR build
+
+# ReMKiT1D Requires gfortran-11, gcc-11 and g++-11. CMake Installs the Package to the 'installs' Folder
+RUN cmake .. -DCMAKE_Fortran_COMPILER=gfortran-11 -DCMAKE_INSTALL_PREFIX=/home/installs/json-fortran -DSKIP_DOC_GEN=TRUE
+RUN make install
+
+WORKDIR /home
+
+# HDF5 is Required and HDF5-1.13.0 is used within the main.yml Workflow Script
+RUN git clone -b hdf5-1_13_0 https://github.com/HDFGroup/hdf5.git
+
+WORKDIR hdf5
+RUN mkdir build
+WORKDIR build
+
+# ReMKiT1D Requires gfortran-11, gcc-11 and g++-11. CMake Installs the Package to the 'installs' Folder
+RUN cmake .. -DCMAKE_Fortran_COMPILER=gfortran-11 -DCMAKE_C_COMPILER=gcc-11 -DCMAKE_CXX_COMPILER=g++-11 -DHDF5_BUILD_FORTRAN=ON -DCMAKE_INSTALL_PREFIX=/home/installs/hdf5
+
+RUN make install
+WORKDIR /home
+
+# ReMKiT1D v1.2.0+ requires sundials 
+RUN git clone -b v7.2.1 https://github.com/LLNL/sundials.git
+WORKDIR sundials
+RUN mkdir build
+WORKDIR build
+
+RUN cmake .. -DENABLE_MPI=ON -DCMAKE_INSTALL_PREFIX=/home/installs/sundials -DENABLE_LAPACK=ON -DBUILD_FORTRAN_MODULE_INTERFACE=ON -DCMAKE_C_COMPILER=gcc-11
+RUN make
+RUN make install
+
+WORKDIR /home
+
+RUN echo LD_LIBRARY_PATH=/home/mpich-install/lib:$LD_LIBRARY_PATH
+
+# Add the PATH to assist CMake and Make Find the Required Compilation Files
+ENV PFUNIT_DIR=/home/installs/pFUnit
+ENV PETSC_DIR=/home/petsc
+ENV PETSC_ARCH=arch-linux-c-opt
+ENV PATH=$PATH:/home/installs/petsc
+ENV PATH=$PATH:/home/installs/hdf5
+ENV PATH=$PATH:/home/installs/json-fortran/jsonfortran-gnu-8.2.5
+ENV PATH=$PATH:/home/mpich-install/bin:$PATH
+ENV PATH=$PATH:/home/installs/sundials
+ENV LD_LIBRARY_PATH=/home/mpich-install/lib:$LD_LIBRARY_PATH
+ENV LD_LIBRARY_PATH=/home/sundials/lib:$LD_LIBRARY_PATH
+ENV SUNDIALS_DIR=/home/installs/sundials/lib/cmake/sundials

docs-pages/01-install.md

@@ -1,6 +1,6 @@
 Title: Building and Installing
 Author: Stefan Mijin 
-Date: 19.09.2024.
+Date: 27.02.2025.
 
 ## Prerequisites
 

docs-pages/02-run.md

@@ -1,6 +1,6 @@
 Title: Running
 Author: Stefan Mijin
-Date: 19.09.2024.
+Date: 25.02.2025.
 
 If the build process is successful, the executable will be in build/src/executables/ReMKiT1D.
 ReMKiT1D requires a config.json file as input. Assuming a valid config file is available, the code can be run from the executable directory using, for example:

docs-pages/index.md

@@ -1,9 +1,9 @@
 Title: Documentation
 Author: Stefan Mijin
-Date: 19.09.2024.
+Date: 27.02.2025.
 
 # **Re**active **M**ultifluid and **Ki**netic **T**ransport in **1D**
 
-rEMKiT1D is a framework for building 1D multi-fluid models of the tokamak Scrape-Off Layer with support for kinetic electron effects and collisional-radiative modelling. The framework is controlled through a Python interface available [here](https://github.com/ukaea/ReMKiT1D-Python), where users can find high-level examples and documentation.
+ReMKiT1D is a framework for building 1D multi-fluid models of the tokamak Scrape-Off Layer with support for kinetic electron effects and collisional-radiative modelling. The framework is controlled through a Python interface available [here](https://github.com/ukaea/ReMKiT1D-Python), where users can find high-level examples and documentation.
 
 For an overview of both the Fortran framework and the Python interface see the [code paper](https://www.sciencedirect.com/science/article/pii/S0010465524001188).

src/modules/abstract_model_builder/term_generator_abstract.f90

@@ -48,6 +48,7 @@ module term_generator_abstract_class
 
         procedure ,public :: getNumImplicitTerms
         procedure ,public :: getNumGeneralTerms
+        procedure ,public :: getGeneratorPrefix
 
         procedure ,public :: moveTerms
 
@@ -110,6 +111,14 @@ pure module function getNumGeneralTerms(this) result(numTerms)
         integer(ik)                         :: numTerms
 
     end function getNumGeneralTerms
+!-----------------------------------------------------------------------------------------------------------------------------------
+    pure module function getGeneratorPrefix(this) result(prefix)
+        !! Get size of this%generalTerms
+
+        class(TermGenerator)   ,intent(in)  :: this
+        character(:) ,allocatable           :: prefix
+
+    end function getGeneratorPrefix
 !-----------------------------------------------------------------------------------------------------------------------------------
     module subroutine moveTerms(this,modelObj,impTermImpGroups,impTermGenGroups,genTermGroups) 
         !! Move terms to modelObj 
@@ -126,4 +135,4 @@ end subroutine moveTerms
 !-----------------------------------------------------------------------------------------------------------------------------------
  end module term_generator_abstract_class
 !-----------------------------------------------------------------------------------------------------------------------------------
-
+

src/modules/abstract_model_builder/term_generator_abstract_procedures.f90

@@ -82,6 +82,16 @@ pure module function getNumGeneralTerms(this) result(numTerms)
     numTerms = size(this%generalTerms)
 
 end function getNumGeneralTerms
+!-----------------------------------------------------------------------------------------------------------------------------------
+    pure module function getGeneratorPrefix(this) result(prefix)
+        !! Get size of this%generalTerms
+
+        class(TermGenerator)   ,intent(in)  :: this
+        character(:) ,allocatable           :: prefix
+
+        prefix = this%generatorPrefix
+
+    end function getGeneratorPrefix
 !-----------------------------------------------------------------------------------------------------------------------------------
 module subroutine moveTerms(this,modelObj,impTermImpGroups,impTermGenGroups,genTermGroups) 
     !! Move terms to modelObj 

src/modules/basic_integrators/implicit_PicardBDE_integrator.f90

@@ -43,11 +43,11 @@ module implicit_PicardBDE_integrator_class
 
         integer(ik) :: maxRestarts = 3 !! Maximum number of consecutive solver restart attempts before critical failure is announced
 
+        integer(ik) :: consolidationInterval = 50 !! How many integration calls before currentNumSubsteps is again reduced to 1
+
         integer(ik) :: hardMaxRestarts = 10 !! Maximum number of consecutive solver restart attempts before critical failure is announced regardless of whether consolidation happened or not
 
         integer(ik) :: restartCount = 0 !! Counter for consecutive number of solver restars 
-
-        integer(ik) :: consolidationInterval = 50 !! How many integration calls before currentNumSubsteps is again reduced to 1
         integer(ik) :: stepsSinceLastConsolidation = 0 !! Counter for steps since last consolidation to 1 substep
 
     end type InternalControllerOptions

src/modules/basic_integrators/implicit_PicardBDE_integrator_procedures.f90

@@ -218,7 +218,6 @@ function checkConvergence(oldVars,newVars,indicesToCheck,nonlinTol,absTol,use2No
 
         end if
         varConverged(i) = relError < nonlinTol .or. absError < epsilon(absError)*absTol
-        !if (.not. varConverged(i)) print*,i,varConverged(i),relError,absError
         if (.not. varConverged(i)) convergenceCounter(i) = convergenceCounter(i) + 1
     end do
 
@@ -354,6 +353,9 @@ subroutine tryIntegrate(this,manipulatedModeller,outputVars,inputVars,numSteps,d
 
         allocate(implicitVectorInit,source=this%implicitVectorOld)
         do i = 1, numSteps
+
+            call printNamedValue(this%integratorName//": Current number of substeps",numSteps)
+            call printNamedValue(this%integratorName//": Starting substep",i)
             if (inputVars%isVarNameRegistered("time")) &
             this%buffer%variables(timeVarIndex)%entry(1) = this%buffer%variables(timeVarIndex)%entry(1) + dt(i)
             tolReached = .false.

src/modules/basic_support/physics_functions.f90

@@ -33,13 +33,14 @@ pure module function logLee(Te,ne) result(res)
 
     end function logLee
 !-----------------------------------------------------------------------------------------------------------------------------------
-    pure module function logLei(Te,ne,Z) result(res)
+    pure module function logLei(Te,ne,Z,removeDisc) result(res)
         !! Calculate Coulomb logarithm for electron-ion collisions (NRL Formulary 2013 page 34 equation b). Assumes  Te > Ti*Z*me/mi 
 
-        real(rk) ,intent(in) :: Te !! Electron temperature in eV
-        real(rk) ,intent(in) :: ne !! Electron density in m^{-3} 
-        real(rk) ,intent(in) :: Z !! Ion charge
-        real(rk)             :: res
+        real(rk) ,intent(in)           :: Te !! Electron temperature in eV
+        real(rk) ,intent(in)           :: ne !! Electron density in m^{-3} 
+        real(rk) ,intent(in)           :: Z !! Ion charge
+        logical  ,intent(in), optional :: removeDisc !! Remove the discontinuity at 10eV by moving the branch switch to Te = Z**2 * e**2 eV. Defaults to .false.
+        real(rk)                       :: res
 
     end function logLei
 !-----------------------------------------------------------------------------------------------------------------------------------
@@ -135,4 +136,4 @@ end function simpleMoment
     end interface
 !-----------------------------------------------------------------------------------------------------------------------------------
 end module physics_functions
-!-----------------------------------------------------------------------------------------------------------------------------------
+!-----------------------------------------------------------------------------------------------------------------------------------

src/modules/basic_support/physics_functions_procedures.f90

@@ -32,15 +32,23 @@ pure module function logLee(Te,ne) result(res)
 
 end function logLee
 !-----------------------------------------------------------------------------------------------------------------------------------
-pure module function logLei(Te,ne,Z) result(res)
+pure module function logLei(Te,ne,Z,removeDisc) result(res)
     !! Calculate Coulomb logarithm for electron-ion collisions (NRL Formulary 2013 page 34 equation b). Assumes Te > TiZm_e/m_i 
 
-    real(rk) ,intent(in) :: Te !! Electron temperature in eV
-    real(rk) ,intent(in) :: ne !! Electron density in  m^{-3} 
-    real(rk) ,intent(in) :: Z !! Ion charge
+    real(rk) ,intent(in)           :: Te !! Electron temperature in eV
+    real(rk) ,intent(in)           :: ne !! Electron density in m^{-3} 
+    real(rk) ,intent(in)           :: Z !! Ion charge
+    logical  ,intent(in), optional :: removeDisc !! Remove the discontinuity at 10eV by moving the branch switch to Te = Z **2 * e**2 eV. Defaults to .false.
     real(rk)             :: res
 
-    if (Te < 10*Z**2) then 
+    real(rk) :: cutoff
+
+    cutoff = 10*Z**2
+    if (present(removeDisc)) then 
+        if (removeDisc) cutoff = exp(real(2,kind=rk))*Z**2
+    end if
+
+    if (Te < cutoff) then 
 
         res = real(23.0d0 - log(sqrt(ne * 1.00D-6) * Z * Te ** (-3.0d0/2.0d0)),kind=rk)
 

src/modules/basic_support/unary_transforms.f90

@@ -266,7 +266,7 @@ pure module function unaryMinmod(input,realParams,intParams,logicalParams) resul
 
     end function unaryMinmod
 !-----------------------------------------------------------------------------------------------------------------------------------
-    pure module function absFloor(input,realParams,intParams,logicalParams) result(output)
+    pure module function unaryAbsFloor(input,realParams,intParams,logicalParams) result(output)
     !! Floor value filter unary wrapper, uses the absolute value of the first real param as the floor value.
     !! The result it sign(input) * max(abs(input),abs(realParams(1)))
 
@@ -276,7 +276,30 @@ pure module function absFloor(input,realParams,intParams,logicalParams) result(o
         logical      ,optional ,dimension(:) ,intent(in) :: logicalParams
         real(rk) ,allocatable ,dimension(:)              :: output
 
-    end function absFloor
+    end function unaryAbsFloor
+!-----------------------------------------------------------------------------------------------------------------------------------
+    pure module function unaryStep(input,realParams,intParams,logicalParams) result(output)
+    !! Unary step function, takes no params. Returns 1 where input is > 0, else 0.
+
+        real(rk)               ,dimension(:) ,intent(in) :: input 
+        real(rk)     ,optional ,dimension(:) ,intent(in) :: realParams
+        integer(ik)  ,optional ,dimension(:) ,intent(in) :: intParams
+        logical      ,optional ,dimension(:) ,intent(in) :: logicalParams
+        real(rk) ,allocatable ,dimension(:)              :: output
+
+    end function unaryStep
+!-----------------------------------------------------------------------------------------------------------------------------------
+    pure module function unaryFilter(input,realParams,intParams,logicalParams) result(output)
+    !! Unary filter function, takes in two realParams, and returns the input, except where input<realParams(1) or
+    !! input>realParams(2)
+
+        real(rk)               ,dimension(:) ,intent(in) :: input 
+        real(rk)     ,optional ,dimension(:) ,intent(in) :: realParams
+        integer(ik)  ,optional ,dimension(:) ,intent(in) :: intParams
+        logical      ,optional ,dimension(:) ,intent(in) :: logicalParams
+        real(rk) ,allocatable ,dimension(:)              :: output
+
+    end function unaryFilter
 !-----------------------------------------------------------------------------------------------------------------------------------
     end interface
 !-----------------------------------------------------------------------------------------------------------------------------------