Millepede II program, subroutines. More...
Go to the source code of this file.
Functions/Subroutines | |
program | mptwo |
Millepede II main program Pede. | |
subroutine | solglo (ivgbi) |
Error for single global parameter from MINRES. | |
subroutine | addcst |
Add constraint information to matrix and vector. | |
subroutine | feasma |
Matrix for feasible solution. | |
subroutine | feasib (concut, iact) |
Make parameters feasible. | |
subroutine | peread (more) |
Read (block of) records from binary files. | |
subroutine | peprep (mode) |
Prepare records. | |
subroutine | isjajb (nst, is, ja, jb, jsp) |
Decode Millepede record. | |
subroutine | loopn |
Loop with fits and sums. | |
subroutine | ploopa (lunp) |
Print title for iteration. | |
subroutine | ploopb (lunp) |
Print iteration line. | |
subroutine | ploopc (lunp) |
Print sub-iteration line. | |
subroutine | ploopd (lunp) |
Print solution line. | |
subroutine | explfc (lunit) |
Print explanation of iteration table. | |
subroutine | mupdat (i, j, add) |
Update element of global matrix. | |
subroutine | loopbf (nrej, ndfs, sndf, dchi2s, numfil, naccf, chi2f, ndff) |
Loop over records in read buffer (block), fits and sums. | |
subroutine | prtglo |
Print final log file. | |
subroutine | avprod (n, x, b) |
Product symmetric matrix times vector. | |
integer(kind=large) function | ijadd (itema, itemb) |
Index for sparse storage. | |
subroutine | sechms (deltat, nhour, minut, secnd) |
Time conversion. | |
INTEGER function | inone (item) |
Translate labels to indices (for global parameters). | |
subroutine | upone |
Update, redefine hash indices. | |
INTEGER function | iprime (n) |
largest prime number < N. | |
subroutine | loop1 |
First data loop (get global labels). | |
subroutine | loop2 |
Second data loop (number of derivatives, global label pairs). | |
subroutine | vmprep (msize) |
Prepare storage for vectors and matrices. | |
subroutine | minver |
Solution by matrix inversion. | |
subroutine | mdiags |
Solution by diagonalization. | |
subroutine | zdiags |
Covariance matrix for diagonalization. | |
subroutine | mminrs |
Solution with MINRES. | |
subroutine | mcsolv (n, x, y) |
Solution for zero band width preconditioner. | |
subroutine | mvsolv (n, x, y) |
Solution for finite band width preconditioner. | |
subroutine | xloopn |
Standard solution algorithm. | |
subroutine | filetc |
Interprete command line option, steering file. | |
subroutine | filetx |
Interprete text files. | |
INTEGER function | nufile (fname) |
Inquire on file. | |
subroutine | intext (text, nline) |
Interprete text. | |
subroutine | additem (length, list, label, value) |
add item to list | |
subroutine | mstart (text) |
Start of 'module' printout. | |
subroutine | mend |
End of 'module' printout. | |
subroutine | mvopen (lun, fname) |
Open file. | |
subroutine | petime |
Print times. | |
subroutine | addsum (add) |
Accurate summation. | |
subroutine | getsum (asum) |
Get accurate sum. |
Millepede II program, subroutines.
Definition in file pede.f90.
subroutine addcst | ( | ) |
subroutine additem | ( | integer, intent(inout) | length, |
type(listitem), dimension(:), intent(inout), allocatable | list, | ||
integer, intent(in) | label, | ||
real, intent(in) | value | ||
) |
subroutine addsum | ( | double precision | add | ) |
subroutine avprod | ( | integer, intent(in) | n, |
double precision, dimension(n), intent(in) | x, | ||
double precision, dimension(n), intent(out) | b | ||
) |
Product symmetric matrix times vector.
A(sym) * X => B. Used by MINRES method (Is most CPU intensive part). The matrix A is the global matrix in full symmetric or (compressed) sparse storage.
[in] | n | size of matrix |
[in] | x | vector X |
[in] | b | result vector B |
subroutine explfc | ( | integer | lunit | ) |
subroutine feasib | ( | real, intent(in) | concut, |
integer, intent(out) | iact | ||
) |
subroutine feasma | ( | ) |
subroutine filetc | ( | ) |
Interprete command line option, steering file.
Fetch and interprete command line options, if steering file specified, check file existence (calling NUFILE)
If no steering file specified, check default steering file.
Create test files for command line option '-t'.
Read steering file, print some lines, detect names of text and binary files, check file existence, store all file names.
Open all binary files.
Definition at line 5698 of file pede.f90.
References intext(), matint(), mptest(), mptst2(), mstart(), nufile(), ratext(), and rltext().
Referenced by mptwo().
subroutine filetx | ( | ) |
Interprete text files.
Reset flags and read steering and all other text files. Print some lines from each file.
Store parameter values, constraints and measurements.
Check flags METSOL (method of solution) and MATSTO (matrix storage mode). Set default values for flags, which are undefined.
Parameter values, format:
1 label 2 (initial) parameter value 3 pre-sigma 4 label 5 (initial) parameter value 6 pre-sigma 7 label ... ... (number of words is multiple of 3)
Constraint and Wconstrained data, format:
1 0 ! constraint header of four words: 2 right-hand-side ! 0 and -1 ... 3 -1; -2 ! ... indicate ... 4 sigma ! ... header 5 label 6 factor 7 label 8 factor 9 ... ... ... (number of words is multiple of 2, at least 6)
Measured data, format:
1 0 ! constraint header of four words: 2 right-hand-side ! 0 and -1 ... 3 -1 ! ... indicate ... 4 sigma ! ... header 5 label 6 factor 7 label 8 factor 9 ... ... ... (number of words is multiple of 2, at least 6)
Definition at line 6119 of file pede.f90.
References intext(), matint(), mend(), and rltext().
Referenced by mptwo().
subroutine getsum | ( | double precision, intent(out) | asum | ) |
integer(kind=large) function ijadd | ( | integer, intent(in) | itema, |
integer, intent(in) | itemb | ||
) |
INTEGER function inone | ( | integer, intent(in) | item | ) |
Translate labels to indices (for global parameters).
Functions INONE and subroutine UPONE are used to collect items, i.e. labels, and to order and translate them.
In the first phase items are collected and stored by calling IRES=INONE(ITEM)
.
At the first entry the two sub-arrays "a" (globalParLabelIndex) and "b" (globalParHashTable) of length 2N are generated with a start length for N=128 entries. In array "a" two words are reserved for each item: (ITEM, count). The function INONE(ITEM) returns the number of the item. At each entry the argument is compared with the already stored items, new items are stored. Search for entries is done using hash-indices, stored in sub-array "b". The initial hash-index is
j = 1 + mod(ITEM, n_prime) + N
where n_prime is the largest prime number less than N. At each entry the count is increased by one. If N items are stored, the size of the sub-arrays is increased by calling CALL UPONE
.
[in] | item | label |
Definition at line 3647 of file pede.f90.
References iprime(), and upone().
Referenced by addcst(), feasib(), feasma(), loop1(), loop2(), loopn(), mvsolv(), peprep(), and xloopn().
subroutine intext | ( | character (len=*), intent(in) | text, |
integer, intent(in) | nline | ||
) |
INTEGER function iprime | ( | integer, intent(in) | n | ) |
subroutine isjajb | ( | integer, intent(in) | nst, |
integer, intent(inout) | is, | ||
integer, intent(out) | ja, | ||
integer, intent(out) | jb, | ||
integer, intent(out) | jsp | ||
) |
Decode Millepede record.
Get indices JA, JB, IS for next measurement within record:
GLDER(JA)
(INDER(JA+J),GLDER(JA+J),J=1,JB-JA-1)
i.e. JB-JA-1 derivativesGLDER(JB)
(INDER(JB+J),GLDER(JB+J),J=1,IS-JB)
i.e. IST-JB derivativesEnd_of_data is indicated by returned values JA=0 and JB=0 Special data are ignored. At end_of_data the info to the special data is returned: IS = pointer to special data; number of words is NSP=-GLDER(IS)
.
[in] | nst | index of last word of record |
[in,out] | is | index of last global derivative (index of first word of record at the first call) |
[out] | ja | index of measured value (=0 at end), = pointer to local derivatives |
[out] | jb | index of standard deviation (=0 at end), = pointer to global derivatives |
[out] | jsp | index to special data |
subroutine loop1 | ( | ) |
First data loop (get global labels).
Read all data files and add all labels to global labels table, add labels from parameters, constraints and measurements (from text files).
Define variable and fixed global parameters (depending on entries and pre-sigma).
Iterate if records had been skipped due to too small read buffer size.
Definition at line 3804 of file pede.f90.
References hmpdef(), inone(), mend(), mstart(), peprep(), peread(), and upone().
Referenced by mptwo().
subroutine loop2 | ( | ) |
Second data loop (number of derivatives, global label pairs).
Calculate maximum number of local, global derivatives and equations per record.
For sparse storage count index pairs with bit (field) counters to construct sparsity structure (row offsets, (compressed) column lists).
Determine read/write cache splitting from average record values (length, global par. vector/matrix).
Check constraints for rank deficit.
Definition at line 4059 of file pede.f90.
References chindl(), ckbits(), clbits(), feasma(), hmpdef(), inbits(), inone(), isjajb(), mend(), mstart(), ndbits(), peprep(), peread(), sort1k(), spbits(), and vmprep().
Referenced by mptwo().
subroutine loopbf | ( | integer, dimension(0:3), intent(inout) | nrej, |
integer, intent(inout) | ndfs, | ||
double precision, intent(inout) | sndf, | ||
double precision, intent(inout) | dchi2s, | ||
integer, intent(in) | numfil, | ||
integer, dimension(numfil), intent(inout) | naccf, | ||
real, dimension(numfil), intent(inout) | chi2f, | ||
integer, dimension(numfil), intent(inout) | ndff | ||
) |
Loop over records in read buffer (block), fits and sums.
Loop over records in current read buffer block (with multiple threads). Perform local fits (optionally with outlier downweigthing) to calculate Chi2, ndf and r.h.s. 'b' of linear equation system A*x=b. In first iteration(s) fill global matrix A.
For the filling of the global matrix each thread creates from his share of local fits (small) udpdate matrices ( from equations (15), (16)) stored in a write buffer. After all events in the read buffer block have been processed the global matrix is being updated from the matrices in the write buffer in parallel (each row by different thread).
The matrices of the local fits are checked for bordered band structure. For border size b and band width m all elements (i,j) are zero for min(i,j)>b and abs(i-j)>m. For sufficient small (b,m) a solution by root free Cholesky decomposition and forward/backward substitution of the band part is much faster compared to inversion (see broken lines in references). Based on the expected computing cost the faster solution method is selected.
[in,out] | nrej | number of rejected records |
[in,out] | ndfs | sum(ndf) |
[in,out] | sndf | sum(weighted ndf) |
[in,out] | dchi2s | sum(weighted chi2) |
[in] | numfil | number of binary files |
[in,out] | naccf | number of accepted records per binary file |
[in,out] | chi2f | sum(chi2/ndf) per binary file |
[in,out] | ndff | sum(ndf) per binary file |
Definition at line 2308 of file pede.f90.
References chindl(), dbavat(), isjajb(), mupdat(), sort1k(), sqmibb(), and sqminv().
Referenced by loopn().
subroutine loopn | ( | ) |
Loop with fits and sums.
Loop over all binary files. Perform local fits to calculate Chi2, ndf and r.h.s. 'b' of linear equation system A*x=b. In first iteration(s) fill matrix A.
Definition at line 1649 of file pede.f90.
References addsum(), getsum(), gmpdef(), hmpdef(), inone(), loopbf(), mupdat(), peprep(), and peread().
Referenced by xloopn().
subroutine mcsolv | ( | integer, intent(in) | n, |
double precision, dimension(n), intent(in) | x, | ||
double precision, dimension(n), intent(out) | y | ||
) |
subroutine mdiags | ( | ) |
subroutine mend | ( | ) |
subroutine minver | ( | ) |
subroutine mminrs | ( | ) |
Solution with MINRES.
Solve A*x=b by minimizing |A*x-b| iteratively. Parallelized (AVPROD).
Use preconditioner with zero (precon) or finite (equdec) band width.
Definition at line 4968 of file pede.f90.
References avprod(), equdec(), mcsolv(), minres(), mvsolv(), and precon().
Referenced by xloopn().
program mptwo | ( | ) |
subroutine mstart | ( | character (len=*), intent(in) | text | ) |
subroutine mupdat | ( | integer, intent(in) | i, |
integer, intent(in) | j, | ||
double precision, intent(in) | add | ||
) |
subroutine mvopen | ( | integer, intent(in) | lun, |
character (len=*), intent(in) | fname | ||
) |
subroutine mvsolv | ( | integer, intent(in) | n, |
double precision, dimension(n), intent(in) | x, | ||
double precision, dimension(n), intent(out) | y | ||
) |
INTEGER function nufile | ( | character (len=*), intent(inout) | fname | ) |
Inquire on file.
Result = 1 for existing binary file, =2 for existing text file, else =0, < 0 open error.
Text file names are recognized by the filename extension, which should contain 'xt' or 'tx'.
[in,out] | fname | file name, optionaly strip prefix. |
Definition at line 6327 of file pede.f90.
References matint().
Referenced by filetc().
subroutine peprep | ( | integer, intent(in) | mode | ) |
subroutine peread | ( | integer, intent(out) | more | ) |
Read (block of) records from binary files.
Optionally using several threads (each file read by single thread). Records larger than the read buffer (ndimbuf
) are skipped. In case of skipped events in the first loop over all binary files the buffer size is adapted to the maximum record size (and the initial loop (LOOP1
) is repeated). C binary files are handled by readc.c and may be gzipped.
[out] | more | more records to come |
The records consist of parallel integer and real arrays:
real array integer array 1 0.0 error count (this record) 2 RMEAS, measured value 0 JA 3 local derivative index of local derivative 4 local derivative index of local derivative 5 ... 6 SIGMA, error (>0) 0 JB global derivative label of global derivative global derivative label of global derivative IST RMEAS, measured value 0 local derivative index of local derivative local derivative index of local derivative ... SIGMA, error 0 global derivative label of global derivative global derivative label of global derivative ... NR global derivative label of global derivative
Definition at line 1207 of file pede.f90.
subroutine petime | ( | ) |
subroutine ploopa | ( | integer, intent(in) | lunp | ) |
subroutine ploopb | ( | integer, intent(in) | lunp | ) |
subroutine ploopc | ( | integer, intent(in) | lunp | ) |
subroutine ploopd | ( | integer, intent(in) | lunp | ) |
subroutine prtglo | ( | ) |
Print final log file.
For each global parameter:
Definition at line 3191 of file pede.f90.
References mvopen().
Referenced by xloopn().
subroutine sechms | ( | real, intent(in) | deltat, |
integer, intent(out) | nhour, | ||
integer, intent(out) | minut, | ||
real, intent(out) | secnd | ||
) |
subroutine solglo | ( | integer, intent(in) | ivgbi | ) |
Error for single global parameter from MINRES.
Calculate single row 'x_i' from inverse matrix by solving A*x_i=b with b=0 except b_i=1.
[in] | ivgbi | index of variable parameter |
Definition at line 767 of file pede.f90.
References avprod(), ijadd(), mcsolv(), minres(), and mvsolv().
Referenced by xloopn().
subroutine upone | ( | ) |
subroutine vmprep | ( | integer(kind=large), dimension(2), intent(in) | msize | ) |
subroutine xloopn | ( | ) |
Standard solution algorithm.
Iterative solution. In current iteration calculate:
ICALCM = +1 Matrix, gradient, Function value & solution ICALCM = 0 gradient, Function value ICALCM = -1 solution ICALCM = -2 end
Solution is obtained by selected method and improved by line search.
Definition at line 5148 of file pede.f90.
References addcst(), dbdot(), feasib(), inone(), loopn(), mdiags(), mend(), minver(), mminrs(), ploopa(), ploopb(), ploopc(), ploopd(), prtglo(), ptldef(), ptline(), ptlprt(), solglo(), and zdiags().
Referenced by mptwo().