# aderrors - Error analysis of Monte Carlo data with Automatic Differentiation `aderrors` is a fortran implementation of the $`\Gamma`$-method for analysis of Monte Carlo data. It uses Automatic Differentiation to perform **exact linear error propagation**. It preforms the computation of gradients and Hessians of arbitrary functions, which allows a robust and **exact error propagation even in iterative algorithms**. - [Features](#features) - [Examples](#examples) - [Simple analysis of MC data](#simple-analysis-of-mc-data) - [A complete example](#a-complete-example) - [A calculator with uncertainties](#a-calculator-with-uncertainties) - [Installation](#installation) - [Requirements](#requirements) - [Instructions](#instructions) - [Using the library](#using-the-library) - [Full documentation](#full-documentation) - [How to cite](#how-to-cite) ## Features - **Exact** linear error propagation, even in iterative algorithms (i.e. error propagation in fit parameters). - Handles data from **any number of ensembles** (i.e. simulations with different parameters). - Support for **replicas** (i.e. several runs with the same simulation parameters). - Irrelagular MC measurements are handled transparently. - Standalone **portable** implementation without any external dependencies. - **Fast** computation of autocorrelation functions with the [FFT package](http://www.kurims.kyoto-u.ac.jp/~ooura/fft.html) (included in the distribution). - **Exact** determination of gradients and Hessians of arbitrary functions. ## Examples This is just a small collection of examples. Note that the basic data type is `uwreal`, that is able to handle MC histories and data with errors. The error is determined by calling the method `uwerr` on the data type. More examples can be found in the `test` directory of the distribution. They are explained in the documentation `doc/aderrors.pdf`. All these examples use the `module simulator` (available in `test/simulator.f90`) to generate autocorrelated data with autocorrelation function $`\Gamma(t) = \sum_k \lambda_k e^{-|t|/\tau_k}`$. ### Simple analysis of MC data This code performs a simple error analysis of the data with the default parameters ($`S_{\tau}=4`$ for automatic windowing, no tail added to the autocorrelation function). The result of the analysis with the $`\Gamma`$-method is compared with the exact values of the error and $`\tau_{\rm int}`$ returned by the `module simulator`. This example is included in the distribution in the file `test/simple.f90`. ```fortran program simple use ISO_FORTRAN_ENV, Only : error_unit, output_unit use numtypes use aderrors use simulator implicit none integer, parameter :: nd = 20000 type (uwreal) :: x real (kind=DP) :: data_x(nd), err, ti, tau(5), lam(5) ! Fill arrays data_x(:) with autocorrelated ! data from the module simulator. tau = (/1.0_DP, 3.0_DP, 4.0_DP, 5.0_DP, 7.34_DP/) lam = (/1.00_DP, 0.87_DP, 1.23_DP, 0.56_DP, 0.87_DP/) call gen_series(data_x, err, ti, tau, lam, 0.3_DP) ! Load data_x(:) measurements in variable x. Use ! default settings (Stau=4, texp=0, 1 replica) x = data_x ! Perform error analysis (optimal window) call x%uwerr() ! Print results write(*,'(1A,1I6,1A)')'** Measurements: ', nd, ' ** ' write(*,100)' - Gamma-method: ', x%value(), " +/- ", x%error(), '( tauint: ', & x%taui(1), " +/- ", x%dtaui(1), ')' write(*,100)' - Exact: ', 0.3_DP, " +/- ", err, '( tauint: ', ti, " )" 100 FORMAT((2X,1A,1F8.5,1A,1F7.5,5X,1A,1F0.2,1A,1F7.5,1A)) stop end program simple ``` Running this code gives as output ``` ** Measurements: 20000 ** - Gamma-method: 0.29883 +/- 0.04271 ( tauint: 4.13 +/- 0.48344) - Exact: 0.30000 +/- 0.04074 ( tauint: 3.82 ) ``` ### A complete example This example is included in the distribution in the file `test/complete.f90`. The analysis is performed for a complicated non-linear function of the two primary observables. It also uses replica for the MC ensemble labeled 1. ```fortran program complete use ISO_FORTRAN_ENV, Only : error_unit, output_unit use numtypes use constants use aderrors use simulator implicit none integer, parameter :: nd = 5000, nrep=4 type (uwreal) :: x, y, z integer :: iflog, ivrep(nrep)=(/1000,30,3070,900/), i, is, ie real (kind=DP) :: data_x(nd), data_y(nd/2), err, ti, texp real (kind=DP) :: tau(4), & lam_x(4)=(/1.0_DP, 0.70_DP, 0.40_DP, 0.40_DP/), & lam_y(4)=(/2.3_DP, 0.40_DP, 0.20_DP, 0.90_DP/) character (len=200) :: flog='history_z.log' ! Fill arrays data_x(:) with autocorrelated ! data from the module simulator. Use nrep replica tau = (/1.0_DP, 3.0_DP, 12.0_DP, 75.0_DP/) texp = maxval(tau) is = 1 do i = 1, nrep ie = is + ivrep(i) - 1 call gen_series(data_x(is:ie), err, ti, tau, lam_x, 0.3_DP) is = ie + 1 end do ! Fill data_y(:) with different values of tau also using ! module simulator forall (i=1:4) tau(i) = real(2*i,kind=DP) call gen_series(data_y, err, ti, tau, lam_y, 1.3_DP) ! Load data_x(:) measurements in variable x. ! Set replica vector, exponential autocorrelation time ! and ensemble ID. x = data_x call x%set_id(1) call x%set_replica(ivrep) call x%set_texp(texp) ! Load data_y(:) measurements in variable y y = data_y call y%set_id(2) ! Exact, transparent error propagation z = sin(x)/(cos(y) + 1.0_DP) ! Attach tail in ensemble with ID 1 when signal in the ! normalized auto-correlation function equals its error call z%set_dsig(1.0_DP,1) ! Set Stau=3 for automatic window in ensemble with ID 2 call z%set_stau(3.0_DP,2) ! Perform error analysis (tails, optimal window,...) call z%uwerr() ! Print results and output details to flog write(*,'(1A,1F8.5,1A,1F8.5)')'** Observable z: ', z%value(), " +/- ", z%error() do i = 1, z%neid() write(*,'(3X,1A,1I3,3X,1F5.2,"%")',advance="no")'Contribution to error from ensemble ID', & z%eid(i), 100.0_DP*z%error_src(i) write(*,'(2X,1A,1F0.4,1A,1F8.4,1A)')'(tau int: ', z%taui(i), " +/- ", z%dtaui(i), ")" end do open(newunit=iflog, file=trim(flog)) call z%print_hist(iflog) close(iflog) stop end program complete ``` Running this code gives as output ``` ** Observable z: 0.24426 +/- 0.05374 Contribution to error from ensemble ID 1 83.93% (tau int: 5.8333 +/- 2.0772) Contribution to error from ensemble ID 2 16.07% (tau int: 2.5724 +/- 0.5268) ``` The file `history_z.log` contains the MC histories, normalized autocorrelation functions and $`\tau_{\rm int}`$. The format is a simple text file: it can be processed by shell scripts (see `tools/plot/plot_hist.sh`) to produce the following graphics. ID1-hist
ID1-hist ID1-hist ### A calculator with uncertainties The module also handles simple data with errors (i.e. **not** a MC history). This is extremely useful since in many occasions we have to use data from statistically independent sources where only the central value and the error is available. With the module `aderrors` these kind of data in handled transparently, just as if it were another ensemble (see `test/calculator.f90`). ```fortran program calculator use aderrors implicit none type (uwreal) :: x, y, z, t x = (/1.223_8, 0.012_8/) ! x = 1.223 +/- 0.012 y = sin(2.0_8*x) z = 1.0_8 + 2.0_8 * sin(x)*cos(x) t = z - y call t%uwerr() write(*,'(1A,1F18.16,1A,1F18.16)')'Exactly one: ', t%value(), " +/- ", t%error() ! 1.0 +/- 0.0 stop end program calculator ``` Running this code gives as output ``` Exactly one: 1.0000000000000000 +/- 0.0000000000000000 ``` ## Installation ### Requirements The code is strict `fortran 2008` compliant. Any compiler that supports this standard can be used. The code has been tested with `gfortran 6.X`, `gfortran 7.X`, `gfortran 8.X`, `intel fortran 17`, `intel fortran 18`. Note that `gfortran 4.X` and `gfortran 5.X` do not support `submodules` (part of the `fortran 2008` standard). This code will not work whith these versions. ### Instructions 1. Download or clone the repository. 1. Edit the `Makefile` in the `build` directory. Change the compiler command/options (variables `FC` and `FOPT`). 1. Compile the library with `GNU make`. 1. Optionally build/run the test codes with `make test`. Executabes will be placed in the `test` directory. 1. If preferred, move the contents of the `include` and `lib` directories somewhere else. ### Using the library 1. Compile your programs with `-I /include`. 1. Link your programs with the `-L /lib` and `-laderr` options. ## Full documentation Look into the `doc/aderrors.pdf` file. ## How to cite If you use this package for a scientific publication, please cite the original work: "Automatic differentiation for error analysis of Monte Carlo data" Alberto Ramos. [arXiv:1809.01289](https://arxiv.org/abs/1809.01289)