Software repository for the reproduction of the test cases from

**Thomas Richter, Rolf Ulrich, Markus Janczyk:**
*Diffusion models with time-dependent parameters: "An analysis of computational effort and accuracy of different numerical methods"*

This software is used in particular for the reproducibility of the results. 

However, the algorithms can also be used directly for own purposes. If you have any questions about possibly necessary adaptations, please contact thomas.richter@ovgu.de.


# Parts of this repository

## General setup

**Python** collects all Python script. Here, **Python/PythonTools** are several internal functions, e.g. the realizations of KFE and random walks. **Python/results** and **Python/pics** are the directories where the results (figures and text-files) are put.

**C++** collects the C++ scripts.



## Case I

Reproduces Case I of the paper (time-independent)

- Python/TestCase1.py

runs the test-case with random walks, integral equation method and with KFE. It produces output in **Python/pics** and **Python/results**. These results will be used in **C++/testcase1.cc** (as reference solution) and by **Python/TestCase1-Plot.py**

- C++/testcase1.cc

runs the stochastic Euler simulation. Script is started by **C++/run-testcase1.sh**. It reads in the reference solution generated by **Python/TestCase1.py** for computing errors. 
- Python/TestCase1-Plot.py

produces Fig. 6 of the paper. It requires the outputs of **Python/TestCase1.py** and **C++/testcase1.cc**

## Case II

Reproduces Case II of the paper (time-dependent thresholds and drift)

- Python/TestCase2.py

runs the test-case with random walks, integral equation method and with KFE. It produces output in **Python/pics** and **Python/results**. These results will be used in **C++/testcase2.cc** (as reference solution) and by **Python/TestCase2-Plot.py**

- C++/testcase2.cc

runs the stochastic Euler simulation. Script is started by **C++/run-testcase2.sh**. It reads in the reference solution generated by **Python/TestCase2.py** for computing errors. 

- Python/TestCase2-Plot.py

produces Fig. 7 of the paper. It requires the outputs of **Python/TestCase2.py** and **C++/testcase2.cc**

- Python/TestCase2-AdjustRandomWalks.py

runs simulations to reproduce Fig. 11 of the paper and implements the modification of the random walk strategy to limit oscillations.

## Case III

Reproduces Case III of the paper (dependency of the accuracy on the derivative of the drift)

- Python/TestCase3.py

runs the test-case with random walks, integral equation method and with KFE for a fixed discretization but with different values of the drift tau. It produces first part of Fig. 8.

- C++/testcase3.cc

runs the stochastic Euler simulation. Script is started by **C++/run-testcase3.sh**. It reads in the reference solution generated by **Python/TestCase3.py** for computing errors. 


- Python/TestCase3-Plot.py

produces second part of Fig. 8. Depends on the output of **Python/TestCase3.py**


## Case IV

Reproduces Case IV of the paper (accuracy and efficiency for Dirac initial data)

- Python/TestCase4.py

runs the test-case with random walks, integral equation and with KFE for a refined discretizations. 

- Python/TestCase4-Plot.py

produces Fig. 9. Depends on the output of **Python/TestCase4.py**

- Python/TestCase4-showsolution.py

Solves with the KFE and plots the solution as surface plot over time and space variable. This skript is used to create Fig. 10 of the paper. Problem parameters and discretization can be adjusted at the top of the script. To test the different stabilization strategies, one can either adjust the value of theta, or one activates Rannacher time-marching by commenting in the marked lines in the skript PythonTools/kfe.py, here in kfe_ale(..) 

## Data Fitting

Python scripts to fit the KFE model to the Data published by Rolf Ulrich et al. in 

**R. Ulrich, H. Schröter, H. Leuthold, T. Birngruber**
*Automatic and controlled stimulus processing in conflict tasks: Superimposed diffusion processes and delta functions.* Cognitive Psychology, 78 , 148–174

- Python/DataFitting-Simon.py

runs the parameter fitting for the Simon task and produces data for Fig. 9 and Table 1.

- Python/Eriksen-Fletcher.py

runs the parameter fitting for the Eriksen Fletcher task and produces data for Fig. 9 and Table 2.

# Installation & running the examples

## Python

The python skripts can just be started. Just note that they depend on each other, i.e.: **Python/TestCase1.py** produces a reference solution that is required by **C++/testcase1.cc** and the results of both are needed in **Python/TestCase1-Plot.py**

The scripts only depend on standard packages like numpy or scipy and all Python environments should work. One suggestion is to use Spyder as part of Anaconda.

## C++

The C++-programs are not intended for performing the simulations in a stand-alone application. Instead, the SDE is simulated for a given number of trials **N_tr** and a given time step **dt** and this simulation is repeated **64** times in order to estimate the average error. It should however be simple to use the scripts as basis for an efficient parallel simulation tool that uses multithreading.


### Configuration
The C++ test cases must be compiled. The test cases are set up to use **cmake**. We suggest the following (in a Linux-environment or on a Mac using homebrew or MacPorts):


1. Create a directory for compilation, e.g. **C++/bin** now called the **bin-dir**
2. In the **bin-dir** calls cmake by **cmake ..** (adjust the path, if the **bin-dir** is not a subdirectory of the **C++-dir**. 
3. Several options can be adjusted. In **C++/bin** call **ccmake .** to make all necessary changes.

If you change the location of the **bin-dir** you will have to modify the run-scripts **run-testcase[123].sh**. 

### Compilation

Initially and whenever you change the code, the programs must be re-compiled

1. In **C++/bin** just call **make**

### Running the examples

The programs are started in **C++**. For each of the test-case there is a skript to start the program. 

1. In **C++** call **sh ./run-testcase1.sh** (or **sh ./run-testcase2.sh**, etc.)

Each script will start the programs several times. For **Case I**, **Case II** and **Case IV** the simulation is started on a sequence of finer and finer discretizations, for **Case III** the value of *tau* will be changed.

The scripts store the output in **C++/results**. Old outputs will be overwritten! Further, the scripts read information about the reference solution from **Python/resuts**. 

The C++ programs use multithreading the OpenMP. If you do not specify the number of threads to be used, all available threads are taken including all hyperthreads. This is usually not efficient it is therefore advisable to set the number of threads by hand, e.g. by calling

**export OMP_NUM_THREADS=8** 

before calling the run-scripts.


# License Information

Initially the software has been written Thomas Richter, Otto-von-Guericke University Magdeburg, Germany in 2022, 2023 (thomas.richter@ovgu.de)

You are free to use the scripts under the *Creative Commons Attribution 4.0 License*.