lxsd / fvn

%\documentclass[a4paper,10pt]{article}

1

%\documentclass[a4paper,10pt]{article}

\documentclass[a4paper,english]{article}

2

\documentclass[a4paper,english]{article}

3

\usepackage[utf8]{inputenc}

4

\usepackage[utf8]{inputenc}

\usepackage{a4wide}

5

\usepackage{a4wide}

\usepackage{eurosym}

6

\usepackage{eurosym}

\usepackage{url}

7

\usepackage{url}

\usepackage[colorlinks=true,hyperfigures=true]{hyperref}

8

\usepackage[colorlinks=true,hyperfigures=true]{hyperref}

%\usepackage{aeguill}

9

%\usepackage{aeguill}

10

\usepackage{graphicx}

11

\usepackage{graphicx}

\usepackage{babel}

12

\usepackage{babel}

\makeatother

13

\makeatother

14

15

%opening

16

%opening

\title{FVN Documentation}

17

\title{FVN Documentation}

\author{William Daniau}

18

\author{William Daniau}

19

20

\begin{document}

21

\begin{document}

22

\maketitle

23

\maketitle

24

%\begin{abstract}

25

%\begin{abstract}

26

%\end{abstract}

27

%\end{abstract}

\tableofcontents

28

\tableofcontents

29

\section{Whatis fvn,licence,disclaimer etc}

30

\section{Whatis fvn,licence,disclaimer etc}

\subsection{Whatis fvn}

31

\subsection{Whatis fvn}

fvn is a Fortran95 mathematical module. It provides various usefull subroutine covering linear algebra, numerical integration, least square polynomial, spline interpolation, zero finding, special functions etc.

32

fvn is a Fortran95 mathematical module. It provides various usefull subroutine covering linear algebra, numerical integration, least square polynomial, spline interpolation, zero finding, special functions etc.

33

Most of the work is done by interfacing Lapack \url{http://www.netlib.org/lapack} which means that Lapack and Blas \url{http://www.netlib.org/blas} must be available on your system for linking fvn. If you use an AMD microprocessor, the good idea is to use ACML ( AMD Core Math Library \url{http://developer.amd.com/acml.jsp} which contains an optimized Blas/Lapack. Fvn also contains a slightly modified version of Quadpack \url{http://www.netlib.org/quadpack} for performing the numerical integration tasks. Finally the fnlib library \url{http://www.netlib.org/fn} has been added for special functions.

34

Most of the work is done by interfacing Lapack \url{http://www.netlib.org/lapack} which means that Lapack and Blas \url{http://www.netlib.org/blas} must be available on your system for linking fvn. If you use an AMD microprocessor, the good idea is to use ACML ( AMD Core Math Library \url{http://developer.amd.com/acml.jsp} which contains an optimized Blas/Lapack. Fvn also contains a slightly modified version of Quadpack \url{http://www.netlib.org/quadpack} for performing the numerical integration tasks. Finally the fnlib library \url{http://www.netlib.org/fn} has been added for special functions.

35

This module has been initially written for the use of the ``Acoustic and microsonic'' group leaded by Sylvain Ballandras in the Time and Frequency Department of institute Femto-ST \url{http://www.femto-st.fr/}.

36

This module has been initially written for the use of the ``Acoustic and microsonic'' group leaded by Sylvain Ballandras in the Time and Frequency Department of institute Femto-ST \url{http://www.femto-st.fr/}.

37

\subsection{Licence}

38

\subsection{Licence}

The licence of fvn is free. You can do whatever you want with this code as far as you credit the authors.

39

The licence of fvn is free. You can do whatever you want with this code as far as you credit the authors.

40

\subsubsection*{Authors}

41

\subsubsection*{Authors}

As of the day this manuel is written there's only one author of fvn :\newline

42

As of the day this manuel is written there's only one author of fvn :\newline

William Daniau\newline

43

William Daniau\newline

william.daniau@femto-st.fr\newline

44

william.daniau@femto-st.fr\newline

45

\subsection{Disclaimer}

46

\subsection{Disclaimer}

The usual disclaimer applied : This software is provided AS IS in the hope it will be usefull. Use it at your own risks. The authors should not be taken responsible of anything that may result by the use of this software.

47

The usual disclaimer applied : This software is provided AS IS in the hope it will be usefull. Use it at your own risks. The authors should not be taken responsible of anything that may result by the use of this software.

48

\section{Naming scheme and convention}

49

\section{Naming scheme and convention}

The naming scheme of the routines is as follow :

50

The naming scheme of the routines is as follow :

\begin{verbatim}

51

\begin{verbatim}

fvn_*_name()

52

fvn_*_name()

\end{verbatim}

53

\end{verbatim}

where * can be s,d,c or z.

54

where * can be s,d,c or z.

\begin{itemize}

55

\begin{itemize}

\item s is for single precision real (real,real*4,real(4),real(kind=4))

56

\item s is for single precision real (real,real*4,real(4),real(kind=4))

\item d for double precision real (double precision,real*8,real(8),real(kind=8))

57

\item d for double precision real (double precision,real*8,real(8),real(kind=8))

\item c for single precision complex (complex,complex*8,complex(4),complex(kind=4))

58

\item c for single precision complex (complex,complex*8,complex(4),complex(kind=4))

\item z for double precision complex (double complex,complex*16,complex(8),complex(kind=8))

59

\item z for double precision complex (double complex,complex*16,complex(8),complex(kind=8))

\end{itemize}

60

\end{itemize}

In the following description of subroutines parameters, input parameters are followed by (in), output parameters by (out) and parameters which are used as input and modified by the subroutine are followed by (inout).

61

In the following description of subroutines parameters, input parameters are followed by (in), output parameters by (out) and parameters which are used as input and modified by the subroutine are followed by (inout).

62

For each routine, there is a generic interface (simply remove \verb'_*' in the name), so using the specific routine is not mandatory.

63

For each routine, there is a generic interface (simply remove \verb'_*' in the name), so using the specific routine is not mandatory.

64

\section{Linear algebra}

65

\section{Linear algebra}

The linear algebra routines of fvn are an interface to lapack, which make it easier to use.

66

The linear algebra routines of fvn are an interface to lapack, which make it easier to use.

\subsection{Matrix inversion}

67

\subsection{Matrix inversion}

\begin{verbatim}

68

\begin{verbatim}

call fvn_matinv(d,a,inva,status)

69

call fvn_matinv(d,a,inva,status)

\end{verbatim}

70

\end{verbatim}

\begin{itemize}

71

\begin{itemize}

\item d (in) is an integer equal to the matrix rank

72

\item d (in) is an integer equal to the matrix rank

\item a (in) is a real or complex matrix. It will remain untouched.

73

\item a (in) is a real or complex matrix. It will remain untouched.

\item inva (out) is a real or complex matrix which contain the inverse of a at the end of the routine

74

\item inva (out) is a real or complex matrix which contain the inverse of a at the end of the routine

\item status (out) is an optional integer equal to zero if something went wrong

75

\item status (out) is an optional integer equal to zero if something went wrong

\end{itemize}

76

\end{itemize}

77

\subsubsection*{Example}

78

\subsubsection*{Example}

\begin{verbatim}

79

\begin{verbatim}

program inv

80

program inv

use fvn

81

use fvn

implicit none

82

implicit none

83

real(8),dimension(3,3) :: a,inva

84

real(8),dimension(3,3) :: a,inva

85

call random_number(a)

86

call random_number(a)

a=a*100

87

a=a*100

88

call fvn_matinv(3,a,inva)

89

call fvn_matinv(3,a,inva)

write (*,*) a

90

write (*,*) a

write (*,*)

91

write (*,*)

write (*,*) inva

92

write (*,*) inva

write (*,*)

93

write (*,*)

write (*,*) matmul(a,inva)

94

write (*,*) matmul(a,inva)

end program

95

end program

\end{verbatim}

96

\end{verbatim}

97

98

99

\subsection{Matrix determinants}

100

\subsection{Matrix determinants}

\begin{verbatim}

101

\begin{verbatim}

det=fvn_det(d,a,status)

102

det=fvn_det(d,a,status)

\end{verbatim}

103

\end{verbatim}

\begin{itemize}

104

\begin{itemize}

\item d (in) is an integer equal to the matrix rank

105

\item d (in) is an integer equal to the matrix rank

\item a (in) is a real or complex matrix. It will remain untouched.

106

\item a (in) is a real or complex matrix. It will remain untouched.

\item status (out) is an optional integer equal to zero if something went wrong

107

\item status (out) is an optional integer equal to zero if something went wrong

\end{itemize}

108

\end{itemize}

109

\subsubsection*{Example}

110

\subsubsection*{Example}

\begin{verbatim}

111

\begin{verbatim}

program det

112

program det

use fvn

113

use fvn

implicit none

114

implicit none

115

real(8),dimension(3,3) :: a

116

real(8),dimension(3,3) :: a

real(8) :: deta

117

real(8) :: deta

integer :: status

118

integer :: status

119

call random_number(a)

120

call random_number(a)

a=a*100

121

a=a*100

122

deta=fvn_det(3,a,status)

123

deta=fvn_det(3,a,status)

write (*,*) a

124

write (*,*) a

write (*,*)

125

write (*,*)

write (*,*) "Det = ",deta

126

write (*,*) "Det = ",deta

end program

127

end program

128

\end{verbatim}

129

\end{verbatim}

130

131

132

\subsection{Matrix condition}

133

\subsection{Matrix condition}

\begin{verbatim}

134

\begin{verbatim}

call fvn_matcon(d,a,rcond,status)

135

call fvn_matcon(d,a,rcond,status)

\end{verbatim}

136

\end{verbatim}

\begin{itemize}

137

\begin{itemize}

\item d (in) is an integer equal to the matrix rank

138

\item d (in) is an integer equal to the matrix rank

\item a (in) is a real or complex matrix. It will remain untouched.

139

\item a (in) is a real or complex matrix. It will remain untouched.

\item rcond (out) is a real of same kind as matrix a, it will contain the reciprocal condition number of the matrix

140

\item rcond (out) is a real of same kind as matrix a, it will contain the reciprocal condition number of the matrix

\item status (out) is an optional integer equal to zero if something went wrong

141

\item status (out) is an optional integer equal to zero if something went wrong

\end{itemize}

142

\end{itemize}

143

The reciprocal condition number is evaluated using the 1-norm and is define as in equation \ref{rconddef}

144

The reciprocal condition number is evaluated using the 1-norm and is define as in equation \ref{rconddef}

\begin{equation}

145

\begin{equation}

R = \frac{1}{norm(A)*norm(invA)}

146

R = \frac{1}{norm(A)*norm(invA)}

\label{rconddef}

147

\label{rconddef}

\end{equation}

148

\end{equation}

149

The 1-norm itself is defined as the maximum value of the columns absolute values (modulus for complex) sum as in equation \ref{l1norm}

150

The 1-norm itself is defined as the maximum value of the columns absolute values (modulus for complex) sum as in equation \ref{l1norm}

\begin{equation}

151

\begin{equation}

L1 = max_j ( \sum_i{\mid A(i,j)\mid} )

152

L1 = max_j ( \sum_i{\mid A(i,j)\mid} )

\label{l1norm}

153

\label{l1norm}

\end{equation}

154

\end{equation}

155

\subsubsection*{Example}

156

\subsubsection*{Example}

\begin{verbatim}

157

\begin{verbatim}

program cond

158

program cond

use fvn

159

use fvn

implicit none

160

implicit none

161

real(8),dimension(3,3) :: a

162

real(8),dimension(3,3) :: a

real(8) :: rcond

163

real(8) :: rcond

integer :: status

164

integer :: status

165

call random_number(a)

166

call random_number(a)

a=a*100

167

a=a*100

168

call fvn_d_matcon(3,a,rcond,status)

169

call fvn_d_matcon(3,a,rcond,status)

write (*,*) a

170

write (*,*) a

write (*,*)

171

write (*,*)

write (*,*) "Cond = ",rcond

172

write (*,*) "Cond = ",rcond

end program

173

end program

174

\end{verbatim}

175

\end{verbatim}

176

177

\subsection{Eigenvalues/Eigenvectors}

178

\subsection{Eigenvalues/Eigenvectors}

\begin{verbatim}

179

\begin{verbatim}

call fvn_matev(d,a,evala,eveca,status)

180

call fvn_matev(d,a,evala,eveca,status)

\end{verbatim}

181

\end{verbatim}

\begin{itemize}

182

\begin{itemize}

\item d (in) is an integer equal to the matrix rank

183

\item d (in) is an integer equal to the matrix rank

\item a (in) is a real or complex matrix. It will remain untouched.

184

\item a (in) is a real or complex matrix. It will remain untouched.

\item evala (out) is a complex array of same kind as a. It contains the eigenvalues of matrix a

185

\item evala (out) is a complex array of same kind as a. It contains the eigenvalues of matrix a

\item eveca (out) is a complex matrix of same kind as a. Its columns are the eigenvectors of matrix a : eveca(:,j)=jth eigenvector associated with eigenvalue evala(j).

186

\item eveca (out) is a complex matrix of same kind as a. Its columns are the eigenvectors of matrix a : eveca(:,j)=jth eigenvector associated with eigenvalue evala(j).

\item status (out) is an optional integer equal to zero if something went wrong

187

\item status (out) is an optional integer equal to zero if something went wrong

\end{itemize}

188

\end{itemize}

189

\subsubsection*{Example}

190

\subsubsection*{Example}

\begin{verbatim}

191

\begin{verbatim}

program eigen

192

program eigen

use fvn

193

use fvn

implicit none

194

implicit none

195

real(8),dimension(3,3) :: a

196

real(8),dimension(3,3) :: a

complex(8),dimension(3) :: evala

197

complex(8),dimension(3) :: evala

complex(8),dimension(3,3) :: eveca

198

complex(8),dimension(3,3) :: eveca

integer :: status,i,j

199

integer :: status,i,j

200

call random_number(a)

201

call random_number(a)

a=a*100

202

a=a*100

203

call fvn_matev(3,a,evala,eveca,status)

204

call fvn_matev(3,a,evala,eveca,status)

write (*,*) a

205

write (*,*) a

write (*,*)

206

write (*,*)

do i=1,3

207

do i=1,3

write(*,*) "Eigenvalue ",i,evala(i)

208

write(*,*) "Eigenvalue ",i,evala(i)

write(*,*) "Associated Eigenvector :"

209

write(*,*) "Associated Eigenvector :"

do j=1,3

210

do j=1,3

write(*,*) eveca(j,i)

211

write(*,*) eveca(j,i)

end do

212

end do

write(*,*)

213

write(*,*)

end do

214

end do

215

end program

216

end program

217

\end{verbatim}

218

\end{verbatim}

219

220

\subsection{Sparse solving}

221

\subsection{Sparse solving}

By interfacing Tim Davis's SuiteSparse from university of Florida \url{http://www.cise.ufl.edu/research/sparse/SuiteSparse/} which is a reference for this kind of problems, fvn provides simple subroutines for solving linear sparse systems.

222

By interfacing Tim Davis's SuiteSparse from university of Florida \url{http://www.cise.ufl.edu/research/sparse/SuiteSparse/} which is a reference for this kind of problems, fvn provides simple subroutines for solving linear sparse systems.

223

The provided routines solves the equation $Ax=B$ where A is sparse and given in its triplet form.

224

The provided routines solves the equation $Ax=B$ where A is sparse and given in its triplet form.

225

\begin{verbatim}

226

\begin{verbatim}

call fvn_sparse_solve(n,nz,T,Ti,Tj,B,x,status)

227

call fvn_sparse_solve(n,nz,T,Ti,Tj,B,x,status)

\end{verbatim}

228

\end{verbatim}

\begin{itemize}

229

\begin{itemize}

\item For this family of subroutine the two letters (zl,zi,dl,di) of the specific interface name decribe the arguments's type. z is for complex(8), d for real(8), l for integer(8) and i for integer(4)

230

\item For this family of subroutine the two letters (zl,zi,dl,di) of the specific interface name decribe the arguments's type. z is for complex(8), d for real(8), l for integer(8) and i for integer(4)

\item n (in) is an integer equal to the matrix rank

231

\item n (in) is an integer equal to the matrix rank

\item nz (in) is an integer equal to the number of non-zero elements

232

\item nz (in) is an integer equal to the number of non-zero elements

\item T(nz) (in) is a complex/real array containing the non-zero elements

233

\item T(nz) (in) is a complex/real array containing the non-zero elements

\item Ti(nz),Tj(nz) (in) are the indexes of the corresponding element of T in the original matrix.

234

\item Ti(nz),Tj(nz) (in) are the indexes of the corresponding element of T in the original matrix.

\item B(n) (in) is a complex/real array containing the second member of the equation.

235

\item B(n) (in) is a complex/real array containing the second member of the equation.

\item x(n) (out) is a complex/real array containing the solution

236

\item x(n) (out) is a complex/real array containing the solution

\item status (out) is an integer which contain non-zero is something went wrong

237

\item status (out) is an integer which contain non-zero is something went wrong

\end{itemize}

238

\end{itemize}

239

\subsubsection*{Example}

240

\subsubsection*{Example}

\begin{verbatim}

241

\begin{verbatim}

program test_sparse

242

program test_sparse

243

use fvn

244

use fvn

implicit none

245

implicit none

246

integer(8), parameter :: nz=12

247

integer(8), parameter :: nz=12

integer(8), parameter :: n=5

248

integer(8), parameter :: n=5

complex(8),dimension(nz) :: A

249

complex(8),dimension(nz) :: A

integer(8),dimension(nz) :: Ti,Tj

250

integer(8),dimension(nz) :: Ti,Tj

complex(8),dimension(n) :: B,x

251

complex(8),dimension(n) :: B,x

integer(8) :: status

252

integer(8) :: status

253

A = (/ (2.,0.),(3.,0.),(3.,0.),(-1.,0.),(4.,0.),(4.,0.),(-3.,0.),&

254

A = (/ (2.,0.),(3.,0.),(3.,0.),(-1.,0.),(4.,0.),(4.,0.),(-3.,0.),&

(1.,0.),(2.,0.),(2.,0.),(6.,0.),(1.,0.) /)

255

(1.,0.),(2.,0.),(2.,0.),(6.,0.),(1.,0.) /)

B = (/ (8.,0.), (45.,0.), (-3.,0.), (3.,0.), (19.,0.)/)

256

B = (/ (8.,0.), (45.,0.), (-3.,0.), (3.,0.), (19.,0.)/)

Ti = (/ 1,2,1,3,5,2,3,4,5,3,2,5 /)

257

Ti = (/ 1,2,1,3,5,2,3,4,5,3,2,5 /)

Tj = (/ 1,1,2,2,2,3,3,3,3,4,5,5 /)

258

Tj = (/ 1,1,2,2,2,3,3,3,3,4,5,5 /)

259

!specific routine that will be used here

260

!specific routine that will be used here

!call fvn_zl_sparse_solve(n,nz,A,Ti,Tj,B,x,status)

261

!call fvn_zl_sparse_solve(n,nz,A,Ti,Tj,B,x,status)

call fvn_sparse_solve(n,nz,A,Ti,Tj,B,x,status)

262

call fvn_sparse_solve(n,nz,A,Ti,Tj,B,x,status)

write(*,*) x

263

write(*,*) x

264

end program

265

end program

266

267

program test_sparse

268

program test_sparse

269

use fvn

270

use fvn

implicit none

271

implicit none

272

integer(4), parameter :: nz=12

273

integer(4), parameter :: nz=12

integer(4), parameter :: n=5

274

integer(4), parameter :: n=5

real(8),dimension(nz) :: A

275

real(8),dimension(nz) :: A

integer(4),dimension(nz) :: Ti,Tj

276

integer(4),dimension(nz) :: Ti,Tj

real(8),dimension(n) :: B,x

277

real(8),dimension(n) :: B,x

integer(4) :: status

278

integer(4) :: status

279

A = (/ 2.,3.,3.,-1.,4.,4.,-3.,1.,2.,2.,6.,1. /)

280

A = (/ 2.,3.,3.,-1.,4.,4.,-3.,1.,2.,2.,6.,1. /)

B = (/ 8., 45., -3., 3., 19./)

281

B = (/ 8., 45., -3., 3., 19./)

Ti = (/ 1,2,1,3,5,2,3,4,5,3,2,5 /)

282

Ti = (/ 1,2,1,3,5,2,3,4,5,3,2,5 /)

Tj = (/ 1,1,2,2,2,3,3,3,3,4,5,5 /)

283

Tj = (/ 1,1,2,2,2,3,3,3,3,4,5,5 /)

284

!specific routine that will be used here

285

!specific routine that will be used here

!call fvn_di_sparse_solve(n,nz,A,Ti,Tj,B,x,status)

286

!call fvn_di_sparse_solve(n,nz,A,Ti,Tj,B,x,status)

call fvn_sparse_solve(n,nz,A,Ti,Tj,B,x,status)

287

call fvn_sparse_solve(n,nz,A,Ti,Tj,B,x,status)

write(*,*) x

288

write(*,*) x

289

end program

290

end program

291

292

293

\end{verbatim}

294

\end{verbatim}

295

\subsection{Identity matrix}

296

\subsection{Identity matrix}

\begin{verbatim}

297

\begin{verbatim}

I=fvn_*_ident(n) (*=s,d,c,z)

298

I=fvn_*_ident(n) (*=s,d,c,z)

\end{verbatim}

299

\end{verbatim}

\begin{itemize}

300

\begin{itemize}

\item n (in) is an integer equal to the matrix rank

301

\item n (in) is an integer equal to the matrix rank

\end{itemize}

302

\end{itemize}

This function return the identity matrix of rank n, in the type of the left hand side. No generic interface for this one.

303

This function return the identity matrix of rank n, in the type of the left hand side. No generic interface for this one.

304

305

306

\section{Interpolation}

307

\section{Interpolation}

308

\subsection{Quadratic Interpolation}

309

\subsection{Quadratic Interpolation}

fvn provide function for interpolating values of a tabulated function of 1, 2 or 3 variables, for both single and double precision.

310

fvn provide function for interpolating values of a tabulated function of 1, 2 or 3 variables, for both single and double precision.

\subsubsection{One variable function}

311

\subsubsection{One variable function}

\begin{verbatim}

312

\begin{verbatim}

value=fvn_quad_interpol(x,n,xdata,ydata)

313

value=fvn_quad_interpol(x,n,xdata,ydata)

\end{verbatim}

314

\end{verbatim}

\begin{itemize}

315

\begin{itemize}

\item x is the real where we want to evaluate the function

316

\item x is the real where we want to evaluate the function

\item n is the number of tabulated values

317

\item n is the number of tabulated values

\item xdata(n) contains the tabulated coordinates

318

\item xdata(n) contains the tabulated coordinates

\item ydata(n) contains the tabulated function values ydata(i)=y(xdata(i))

319

\item ydata(n) contains the tabulated function values ydata(i)=y(xdata(i))

\end{itemize}

320

\end{itemize}

xdata must be strictly increasingly ordered.

321

xdata must be strictly increasingly ordered.

x must be within the range of xdata to actually perform an interpolation, otherwise the resulting value is an extrapolation

322

x must be within the range of xdata to actually perform an interpolation, otherwise the resulting value is an extrapolation

\paragraph*{Example}

323

\paragraph*{Example}

\begin{verbatim}

324

\begin{verbatim}

program inter1d

325

program inter1d

326

use fvn

327

use fvn

implicit none

328

implicit none

329

integer(kind=4),parameter :: ndata=33

330

integer(kind=4),parameter :: ndata=33

integer(kind=4) :: i,nout

331

integer(kind=4) :: i,nout

real(kind=8) :: f,fdata(ndata),h,pi,q,sin,x,xdata(ndata)

332

real(kind=8) :: f,fdata(ndata),h,pi,q,sin,x,xdata(ndata)

real(kind=8) ::tv

333

real(kind=8) ::tv

334

intrinsic sin

335

intrinsic sin

336

f(x)=sin(x)

337

f(x)=sin(x)

338

xdata(1)=0.

339

xdata(1)=0.

fdata(1)=f(xdata(1))

340

fdata(1)=f(xdata(1))

h=1./32.

341

h=1./32.

do i=2,ndata

342

do i=2,ndata

xdata(i)=xdata(i-1)+h

343

xdata(i)=xdata(i-1)+h

fdata(i)=f(xdata(i))

344

fdata(i)=f(xdata(i))

end do

345

end do

call random_seed()

346

call random_seed()

call random_number(x)

347

call random_number(x)

348

q=fvn_d_quad_interpol(x,ndata,xdata,fdata)

349

q=fvn_d_quad_interpol(x,ndata,xdata,fdata)

350

tv=f(x)

351

tv=f(x)

write(*,*) "x ",x

352

write(*,*) "x ",x

write(*,*) "Calculated (real) value :",tv

353

write(*,*) "Calculated (real) value :",tv

write(*,*) "fvn interpolation :",q

354

write(*,*) "fvn interpolation :",q

write(*,*) "Relative fvn error :",abs((q-tv)/tv)

355

write(*,*) "Relative fvn error :",abs((q-tv)/tv)

356

end program

357

end program

358

\end{verbatim}

359

\end{verbatim}

360

361

\subsubsection{Two variables function}

362

\subsubsection{Two variables function}

\begin{verbatim}

363

\begin{verbatim}

value=fvn_quad_2d_interpol(x,y,nx,xdata,ny,ydata,zdata)

364

value=fvn_quad_2d_interpol(x,y,nx,xdata,ny,ydata,zdata)

\end{verbatim}

365

\end{verbatim}

\begin{itemize}

366

\begin{itemize}

\item x,y are the real coordinates where we want to evaluate the function

367

\item x,y are the real coordinates where we want to evaluate the function

\item nx is the number of tabulated values along x axis

368

\item nx is the number of tabulated values along x axis

\item xdata(nx) contains the tabulated x

369

\item xdata(nx) contains the tabulated x

\item ny is the number of tabulated values along y axis

370

\item ny is the number of tabulated values along y axis

\item ydata(ny) contains the tabulated y

371

\item ydata(ny) contains the tabulated y

\item zdata(nx,ny) contains the tabulated function values zdata(i,j)=z(xdata(i),ydata(j))

372

\item zdata(nx,ny) contains the tabulated function values zdata(i,j)=z(xdata(i),ydata(j))

\end{itemize}

373

\end{itemize}

xdata and ydata must be strictly increasingly ordered.

374

xdata and ydata must be strictly increasingly ordered.

(x,y) must be within the range of xdata and ydata to actually perform an interpolation, otherwise the resulting value is an extrapolation

375

(x,y) must be within the range of xdata and ydata to actually perform an interpolation, otherwise the resulting value is an extrapolation

376

\paragraph*{Example}

377

\paragraph*{Example}

378

\begin{verbatim}

379

\begin{verbatim}

program inter2d

380

program inter2d

use fvn

381

use fvn

implicit none

382

implicit none

383

integer(kind=4),parameter :: nx=21,ny=42

384

integer(kind=4),parameter :: nx=21,ny=42

integer(kind=4) :: i,j

385

integer(kind=4) :: i,j

real(kind=8) :: f,fdata(nx,ny),dble,pi,q,sin,x,xdata(nx),y,ydata(ny)

386

real(kind=8) :: f,fdata(nx,ny),dble,pi,q,sin,x,xdata(nx),y,ydata(ny)

real(kind=8) :: tv

387

real(kind=8) :: tv

388

intrinsic dble,sin

389

intrinsic dble,sin

390

f(x,y)=sin(x+2.*y)

391

f(x,y)=sin(x+2.*y)

do i=1,nx

392

do i=1,nx

xdata(i)=dble(i-1)/dble(nx-1)

393

xdata(i)=dble(i-1)/dble(nx-1)

end do

394

end do

do i=1,ny

395

do i=1,ny

ydata(i)=dble(i-1)/dble(ny-1)

396

ydata(i)=dble(i-1)/dble(ny-1)

end do

397

end do

do i=1,nx

398

do i=1,nx

do j=1,ny

399

do j=1,ny

fdata(i,j)=f(xdata(i),ydata(j))

400

fdata(i,j)=f(xdata(i),ydata(j))

end do

401

end do

402

end do

call random_seed()

403

call random_seed()

call random_number(x)

404

call random_number(x)

call random_number(y)

405

call random_number(y)

406

q=fvn_d_quad_2d_interpol(x,y,nx,xdata,ny,ydata,fdata)

407

q=fvn_d_quad_2d_interpol(x,y,nx,xdata,ny,ydata,fdata)

tv=f(x,y)

408

tv=f(x,y)

409

write(*,*) "x y",x,y

410

write(*,*) "x y",x,y

write(*,*) "Calculated (real) value :",tv

411

write(*,*) "Calculated (real) value :",tv

write(*,*) "fvn interpolation :",q

412

write(*,*) "fvn interpolation :",q

write(*,*) "Relative fvn error :",abs((q-tv)/tv)

413

write(*,*) "Relative fvn error :",abs((q-tv)/tv)

414

end program

415

end program

416

\end{verbatim}

417

\end{verbatim}

418

419

420

\subsubsection{Three variables function}

421

\subsubsection{Three variables function}

\begin{verbatim}

422

\begin{verbatim}

value=fvn_quad_3d_interpol(x,y,z,nx,xdata,ny,ydata,nz,zdata,tdata)

423

value=fvn_quad_3d_interpol(x,y,z,nx,xdata,ny,ydata,nz,zdata,tdata)

\end{verbatim}

424

\end{verbatim}

\begin{itemize}

425

\begin{itemize}

\item x,y,z are the real coordinates where we want to evaluate the function

426

\item x,y,z are the real coordinates where we want to evaluate the function

\item nx is the number of tabulated values along x axis

427

\item nx is the number of tabulated values along x axis

\item xdata(nx) contains the tabulated x

428

\item xdata(nx) contains the tabulated x

\item ny is the number of tabulated values along y axis

429

\item ny is the number of tabulated values along y axis

\item ydata(ny) contains the tabulated y

430

\item ydata(ny) contains the tabulated y

\item nz is the number of tabulated values along z axis

431

\item nz is the number of tabulated values along z axis

\item zdata(ny) contains the tabulated z

432

\item zdata(ny) contains the tabulated z

\item tdata(nx,ny,nz) contains the tabulated function values tdata(i,j,k)=t(xdata(i),ydata(j),zdata(k))

433

\item tdata(nx,ny,nz) contains the tabulated function values tdata(i,j,k)=t(xdata(i),ydata(j),zdata(k))

\end{itemize}

434

\end{itemize}

xdata, ydata and zdata must be strictly increasingly ordered.

435

xdata, ydata and zdata must be strictly increasingly ordered.

(x,y,z) must be within the range of xdata and ydata to actually perform an interpolation, otherwise the resulting value is an extrapolation

436

(x,y,z) must be within the range of xdata and ydata to actually perform an interpolation, otherwise the resulting value is an extrapolation

437

\paragraph*{Example}

438

\paragraph*{Example}

\begin{verbatim}

439

\begin{verbatim}

program inter3d

440

program inter3d

use fvn

441

use fvn

442

implicit none

443

implicit none

444

integer(kind=4),parameter :: nx=21,ny=42,nz=18

445

integer(kind=4),parameter :: nx=21,ny=42,nz=18

integer(kind=4) :: i,j,k

446

integer(kind=4) :: i,j,k

real(kind=8) :: f,fdata(nx,ny,nz),dble,pi,q,sin,x,xdata(nx),y,ydata(ny),z,zdata(nz)

447

real(kind=8) :: f,fdata(nx,ny,nz),dble,pi,q,sin,x,xdata(nx),y,ydata(ny),z,zdata(nz)

real(kind=8) :: tv

448

real(kind=8) :: tv

449

intrinsic dble,sin

450

intrinsic dble,sin

451

f(x,y,z)=sin(x+2.*y+3.*z)

452

f(x,y,z)=sin(x+2.*y+3.*z)

do i=1,nx

453

do i=1,nx

xdata(i)=2.*(dble(i-1)/dble(nx-1))

454

xdata(i)=2.*(dble(i-1)/dble(nx-1))

end do

455

end do

do i=1,ny

456

do i=1,ny

ydata(i)=2.*(dble(i-1)/dble(ny-1))

457

ydata(i)=2.*(dble(i-1)/dble(ny-1))

end do

458

end do

do i=1,nz

459

do i=1,nz

zdata(i)=2.*(dble(i-1)/dble(nz-1))

460

zdata(i)=2.*(dble(i-1)/dble(nz-1))

end do

461

end do

do i=1,nx

462

do i=1,nx

do j=1,ny

463

do j=1,ny

do k=1,nz

464

do k=1,nz

fdata(i,j,k)=f(xdata(i),ydata(j),zdata(k))

465

fdata(i,j,k)=f(xdata(i),ydata(j),zdata(k))

end do

466

end do

467

end do

468

end do

call random_seed()

469

call random_seed()

call random_number(x)

470

call random_number(x)

call random_number(y)

471

call random_number(y)

call random_number(z)

472

call random_number(z)

473

q=fvn_d_quad_3d_interpol(x,y,z,nx,xdata,ny,ydata,nz,zdata,fdata)

474

q=fvn_d_quad_3d_interpol(x,y,z,nx,xdata,ny,ydata,nz,zdata,fdata)

tv=f(x,y,z)

475

tv=f(x,y,z)

476

write(*,*) "x y z",x,y,z

477

write(*,*) "x y z",x,y,z

write(*,*) "Calculated (real) value :",tv

478

write(*,*) "Calculated (real) value :",tv

write(*,*) "fvn interpolation :",q

479

write(*,*) "fvn interpolation :",q

write(*,*) "Relative fvn error :",abs((q-tv)/tv)

480

write(*,*) "Relative fvn error :",abs((q-tv)/tv)

481

end program

482

end program

483

\end{verbatim}

484

\end{verbatim}

485

\subsubsection{Utility procedure}

486

\subsubsection{Utility procedure}

fvn provides a simple utility procedure to locate the interval in which a value is located in an increasingly ordered array.

487

fvn provides a simple utility procedure to locate the interval in which a value is located in an increasingly ordered array.

\begin{verbatim}

488

\begin{verbatim}

call fvn_find_interval(x,i,xdata,n)

489

call fvn_find_interval(x,i,xdata,n)

\end{verbatim}

490

\end{verbatim}

\begin{itemize}

491

\begin{itemize}

\item x (in) the real value to locate

492

\item x (in) the real value to locate

\item i (out) the resulting indice

493

\item i (out) the resulting indice

\item xdata(n) (in) increasingly ordered array

494

\item xdata(n) (in) increasingly ordered array

\item n (in) size of the array

495

\item n (in) size of the array

\end{itemize}

496

\end{itemize}

The resulting integer i is as : $xdata(i) <= x < xdata(i+1)$. If $x < xdata(1)$ then $i=0$ is returned. If $x > xdata(n)$ then $i=n$ is returned. Finally if $x=xdata(n)$ then $i=n-1$ is returned.

497

The resulting integer i is as : $xdata(i) <= x < xdata(i+1)$. If $x < xdata(1)$ then $i=0$ is returned. If $x > xdata(n)$ then $i=n$ is returned. Finally if $x=xdata(n)$ then $i=n-1$ is returned.

498

499

500

\subsection{Akima spline}

501

\subsection{Akima spline}

fvn provides Akima spline interpolation and evaluation for both single and double precision real.

502

fvn provides Akima spline interpolation and evaluation for both single and double precision real.

\subsubsection{Interpolation}

503

\subsubsection{Interpolation}

\begin{verbatim}

504

\begin{verbatim}

call fvn_akima(n,x,y,br,co)

505

call fvn_akima(n,x,y,br,co)

\end{verbatim}

506

\end{verbatim}

\begin{itemize}

507

\begin{itemize}

\item n (in) is an integer equal to the number of points

508

\item n (in) is an integer equal to the number of points

\item x(n) (in) ,y(n) (in) are the known couples of coordinates

509

\item x(n) (in) ,y(n) (in) are the known couples of coordinates

\item br (out) on output contains a copy of x

510

\item br (out) on output contains a copy of x

\item co(4,n) (out) is a real matrix containing the 4 coefficients of the Akima interpolation spline for a given interval.

511

\item co(4,n) (out) is a real matrix containing the 4 coefficients of the Akima interpolation spline for a given interval.

\end{itemize}

512

\end{itemize}

513

\subsubsection{Evaluation}

514

\subsubsection{Evaluation}

\begin{verbatim}

515

\begin{verbatim}

y=fvn_spline_eval(x,n,br,co)

516

y=fvn_spline_eval(x,n,br,co)

\end{verbatim}

517

\end{verbatim}

\begin{itemize}

518

\begin{itemize}

\item x (in) is the point where we want to evaluate

519

\item x (in) is the point where we want to evaluate

\item n (in) is the number of known points and br(n) (in), co(4,n) (in) \\

520

\item n (in) is the number of known points and br(n) (in), co(4,n) (in) \\

are the outputs of fvn\_x\_akima(n,x,y,br,co)

521

are the outputs of fvn\_x\_akima(n,x,y,br,co)

\end{itemize}

522

\end{itemize}

523

\subsubsection{Example}

524

\subsubsection{Example}

In the following example we will use Akima splines to interpolate a sinus function with 30 points between -10 and 10. We then use the evaluation function to calculate the coordinates of 1000 points between -11 and 11, and write a 3 columns file containing : x, calculated sin(x), interpolation evaluation of sin(x).

525

In the following example we will use Akima splines to interpolate a sinus function with 30 points between -10 and 10. We then use the evaluation function to calculate the coordinates of 1000 points between -11 and 11, and write a 3 columns file containing : x, calculated sin(x), interpolation evaluation of sin(x).

526

One can see that the interpolation is very efficient even with only 30 points. Of course as soon as we leave the -10 to 10 interval, the values are extrapolated and thus can lead to very inacurrate values.

527

One can see that the interpolation is very efficient even with only 30 points. Of course as soon as we leave the -10 to 10 interval, the values are extrapolated and thus can lead to very inacurrate values.

528

\begin{verbatim}

529

\begin{verbatim}

program akima

530

program akima

use fvn

531

use fvn

implicit none

532

implicit none

533

integer :: nbpoints,nppoints,i

534

integer :: nbpoints,nppoints,i

real(8),dimension(:),allocatable :: x_d,y_d,breakpoints_d

535

real(8),dimension(:),allocatable :: x_d,y_d,breakpoints_d

real(8),dimension(:,:),allocatable :: coeff_fvn_d

536

real(8),dimension(:,:),allocatable :: coeff_fvn_d

real(8) :: xstep_d,xp_d,ty_d,fvn_y_d

537

real(8) :: xstep_d,xp_d,ty_d,fvn_y_d

538

open(2,file='fvn_akima_double.dat')

539

open(2,file='fvn_akima_double.dat')

open(3,file='fvn_akima_breakpoints_double.dat')

540

open(3,file='fvn_akima_breakpoints_double.dat')

nbpoints=30

541

nbpoints=30

allocate(x_d(nbpoints))

542

allocate(x_d(nbpoints))

allocate(y_d(nbpoints))

543

allocate(y_d(nbpoints))

allocate(breakpoints_d(nbpoints))

544

allocate(breakpoints_d(nbpoints))

allocate(coeff_fvn_d(4,nbpoints))

545

allocate(coeff_fvn_d(4,nbpoints))

546

xstep_d=20./dfloat(nbpoints)

547

xstep_d=20./dfloat(nbpoints)

do i=1,nbpoints

548

do i=1,nbpoints

x_d(i)=-10.+dfloat(i)*xstep_d

549

x_d(i)=-10.+dfloat(i)*xstep_d

y_d(i)=dsin(x_d(i))

550

y_d(i)=dsin(x_d(i))

write(3,44) (x_d(i),y_d(i))

551

write(3,44) (x_d(i),y_d(i))

end do

552

end do

close(3)

553

close(3)

554

call fvn_d_akima(nbpoints,x_d,y_d,breakpoints_d,coeff_fvn_d)

555

call fvn_d_akima(nbpoints,x_d,y_d,breakpoints_d,coeff_fvn_d)

556

nppoints=1000

557

nppoints=1000

xstep_d=22./dfloat(nppoints)

558

xstep_d=22./dfloat(nppoints)

do i=1,nppoints

559

do i=1,nppoints

xp_d=-11.+dfloat(i)*xstep_d

560

xp_d=-11.+dfloat(i)*xstep_d

ty_d=dsin(xp_d)

561

ty_d=dsin(xp_d)

fvn_y_d=fvn_d_spline_eval(xp_d,nbpoints-1,breakpoints_d,coeff_fvn_d)

562

fvn_y_d=fvn_d_spline_eval(xp_d,nbpoints-1,breakpoints_d,coeff_fvn_d)

write(2,44) (xp_d,ty_d,fvn_y_d)

563

write(2,44) (xp_d,ty_d,fvn_y_d)

end do

564

end do

565

close(2)

566

close(2)

567

44 FORMAT(4(1X,1PE22.14))

568

44 FORMAT(4(1X,1PE22.14))

569

end program

570

end program

571

\end{verbatim}

572

\end{verbatim}

Results are plotted on figure \ref{akima}

573

Results are plotted on figure \ref{akima}

574

\begin{figure}

575

\begin{figure}

\begin{center}

576

\begin{center}

\includegraphics[width=0.9\textwidth]{akima.pdf}

577

\includegraphics[width=0.9\textwidth]{akima.pdf}

% akima.pdf: 504x720 pixel, 72dpi, 17.78x25.40 cm, bb=0 0 504 720

578

% akima.pdf: 504x720 pixel, 72dpi, 17.78x25.40 cm, bb=0 0 504 720

\caption{Akima Spline Interpolation}

579

\caption{Akima Spline Interpolation}

\label{akima}

580

\label{akima}

\end{center}

581

\end{center}

582

\end{figure}

583

\end{figure}

584

585

586

\section{Least square polynomial}

587

\section{Least square polynomial}

fvn provide a function to find a least square polynomial of a given degree, for real in single or double precision. It is performed using Lapack subroutine sgelss (dgelss), which solve this problem using singular value decomposition.

588

fvn provide a function to find a least square polynomial of a given degree, for real in single or double precision. It is performed using Lapack subroutine sgelss (dgelss), which solve this problem using singular value decomposition.

589

\begin{verbatim}

590

\begin{verbatim}

call fvn_lspoly(np,x,y,deg,coeff,status)

591

call fvn_lspoly(np,x,y,deg,coeff,status)

\end{verbatim}

592

\end{verbatim}

\begin{itemize}

593

\begin{itemize}

\item np (in) is an integer equal to the number of points

594

\item np (in) is an integer equal to the number of points

\item x(np) (in),y(np) (in) are the known coordinates

595

\item x(np) (in),y(np) (in) are the known coordinates

\item deg (in) is an integer equal to the degree of the desired polynomial, it must be lower than np.

596

\item deg (in) is an integer equal to the degree of the desired polynomial, it must be lower than np.

\item coeff(deg+1) (out) on output contains the polynomial coefficients

597

\item coeff(deg+1) (out) on output contains the polynomial coefficients

\item status (out) is an integer containing 0 if a problem occured.

598

\item status (out) is an integer containing 0 if a problem occured.

\end{itemize}

599

\end{itemize}

600

\subsection*{Example}

601

\subsection*{Example}

Here's a simple example : we've got 13 measurement points and we want to find the least square degree 3 polynomial for these points :

602

Here's a simple example : we've got 13 measurement points and we want to find the least square degree 3 polynomial for these points :

\begin{verbatim}

603

\begin{verbatim}

program lsp

604

program lsp

use fvn

605

use fvn

implicit none

606

implicit none

607

integer,parameter :: npoints=13,deg=3

608

integer,parameter :: npoints=13,deg=3

integer :: status,i

609

integer :: status,i

real(kind=8) :: xm(npoints),ym(npoints),xstep,xc,yc

610

real(kind=8) :: xm(npoints),ym(npoints),xstep,xc,yc

real(kind=8) :: coeff(deg+1)

611

real(kind=8) :: coeff(deg+1)

612

xm = (/ -3.8,-2.7,-2.2,-1.9,-1.1,-0.7,0.5,1.7,2.,2.8,3.2,3.8,4. /)

613

xm = (/ -3.8,-2.7,-2.2,-1.9,-1.1,-0.7,0.5,1.7,2.,2.8,3.2,3.8,4. /)

ym = (/ -3.1,-2.,-0.9,0.8,1.8,0.4,2.1,1.8,3.2,2.8,3.9,5.2,7.5 /)

614

ym = (/ -3.1,-2.,-0.9,0.8,1.8,0.4,2.1,1.8,3.2,2.8,3.9,5.2,7.5 /)

615

open(2,file='fvn_lsp_double_mesure.dat')

616

open(2,file='fvn_lsp_double_mesure.dat')

open(3,file='fvn_lsp_double_poly.dat')

617

open(3,file='fvn_lsp_double_poly.dat')

618

do i=1,npoints

619

do i=1,npoints

write(2,44) xm(i),ym(i)

620

write(2,44) xm(i),ym(i)

end do

621

end do

close(2)

622

close(2)

623

624

call fvn_d_lspoly(npoints,xm,ym,deg,coeff,status)

625

call fvn_d_lspoly(npoints,xm,ym,deg,coeff,status)

626

xstep=(xm(npoints)-xm(1))/1000.

627

xstep=(xm(npoints)-xm(1))/1000.

do i=1,1000

628

do i=1,1000

xc=xm(1)+(i-1)*xstep

629

xc=xm(1)+(i-1)*xstep

yc=poly(xc,coeff)

630

yc=poly(xc,coeff)

write(3,44) xc,yc

631

write(3,44) xc,yc

end do

632

end do

close(3)

633

close(3)

634

44 FORMAT(4(1X,1PE22.14))

635

44 FORMAT(4(1X,1PE22.14))

636

contains

637

contains

function poly(x,coeff)

638

function poly(x,coeff)

implicit none

639

implicit none

real(8) :: x

640

real(8) :: x

real(8) :: coeff(deg+1)

641

real(8) :: coeff(deg+1)

real(8) :: poly

642

real(8) :: poly

integer :: i

643

integer :: i

644

poly=0.

645

poly=0.

646

do i=1,deg+1

647

do i=1,deg+1

poly=poly+coeff(i)*x**(i-1)

648

poly=poly+coeff(i)*x**(i-1)

end do

649

end do

650

end function

651

end function

end program

652

end program

\end{verbatim}

653

\end{verbatim}

The results are plotted on figure \ref{lsp} .

654

The results are plotted on figure \ref{lsp} .

655

\begin{figure}

656

\begin{figure}

\begin{center}

657

\begin{center}

\includegraphics[width=0.9\textwidth]{lsp.pdf}

658

\includegraphics[width=0.9\textwidth]{lsp.pdf}

\caption{Least Square Polynomial}

659

\caption{Least Square Polynomial}

\label{lsp}

660

\label{lsp}

\end{center}

661

\end{center}

\end{figure}

662

\end{figure}

663

664

665

\section{Zero finding}

666

\section{Zero finding}

fvn provide a routine for finding zeros of a complex function using Muller algorithm (only for double complex type). It is based on a version provided on the web by Hans D Mittelmann \url{http://plato.asu.edu/ftp/other\_software/muller.f}.

667

fvn provide a routine for finding zeros of a complex function using Muller algorithm (only for double complex type). It is based on a version provided on the web by Hans D Mittelmann \url{http://plato.asu.edu/ftp/other\_software/muller.f}.

668

\begin{verbatim}

669

\begin{verbatim}

call fvn_muller(f,eps,eps1,kn,nguess,n,x,itmax,infer,ier)

670

call fvn_muller(f,eps,eps1,kn,nguess,n,x,itmax,infer,ier)

\end{verbatim}

671

\end{verbatim}

\begin{itemize}

672

\begin{itemize}

\item f (in) is the complex function (kind=8) for which we search zeros

673

\item f (in) is the complex function (kind=8) for which we search zeros

\item eps (in) is a real(8) corresponding to the first stopping criterion : let fp(z)=f(z)/p where p = (z-z(1))*(z-z(2))*,,,*(z-z(k-1)) and z(1),...,z(k-1) are previously found roots. if ((cdabs(f(z)).le.eps) .and. (cdabs(fp(z)).le.eps)), then z is accepted as a root.

674

\item eps (in) is a real(8) corresponding to the first stopping criterion : let fp(z)=f(z)/p where p = (z-z(1))*(z-z(2))*,,,*(z-z(k-1)) and z(1),...,z(k-1) are previously found roots. if ((cdabs(f(z)).le.eps) .and. (cdabs(fp(z)).le.eps)), then z is accepted as a root.

\item eps1 (in) is a real(8) corresponding to the second stopping criterion : a root is accepted if two successive approximations to a given root agree within eps1. Note that if either or both of the stopping criteria are fulfilled, the root is accepted.

675

\item eps1 (in) is a real(8) corresponding to the second stopping criterion : a root is accepted if two successive approximations to a given root agree within eps1. Note that if either or both of the stopping criteria are fulfilled, the root is accepted.

\item kn (in) is an integer equal to the number of known roots, which must be stored in x(1),...,x(kn), prior to entry in the subroutine.

676

\item kn (in) is an integer equal to the number of known roots, which must be stored in x(1),...,x(kn), prior to entry in the subroutine.

\item nguess (in) is the number of initial guesses provided. These guesses must be stored in x(kn+1),...,x(kn+nguess). nguess must be set equal to zero if no guesses are provided.

677

\item nguess (in) is the number of initial guesses provided. These guesses must be stored in x(kn+1),...,x(kn+nguess). nguess must be set equal to zero if no guesses are provided.

\item n (in) is an integer equal to the number of new roots to be found.

678

\item n (in) is an integer equal to the number of new roots to be found.

\item x (inout) is a complex(8) vector of length kn+n. x(1),...,x(kn) on input must contain any known roots. x(kn+1),..., x(kn+n) on input may, on user option, contain initial guesses for the n new roots which are to be computed. If the user does not provide an initial guess, zero is used. On output, x(kn+1),...,x(kn+n) contain the approximate roots found by the subroutine.

679

\item x (inout) is a complex(8) vector of length kn+n. x(1),...,x(kn) on input must contain any known roots. x(kn+1),..., x(kn+n) on input may, on user option, contain initial guesses for the n new roots which are to be computed. If the user does not provide an initial guess, zero is used. On output, x(kn+1),...,x(kn+n) contain the approximate roots found by the subroutine.

\item itmax (in) is an integer equal to the maximum allowable number of iterations per root.

680

\item itmax (in) is an integer equal to the maximum allowable number of iterations per root.

\item infer (out) is an integer vector of size kn+n. On output infer(j) contains the number of iterations used in finding the j-th root when convergence was achieved. If convergence was not obtained in itmax iterations, infer(j) will be greater than itmax

681

\item infer (out) is an integer vector of size kn+n. On output infer(j) contains the number of iterations used in finding the j-th root when convergence was achieved. If convergence was not obtained in itmax iterations, infer(j) will be greater than itmax

\item ier (out) is an integer used as an error parameter. ier = 33 indicates failure to converge within itmax iterations for at least one of the (n) new roots.

682

\item ier (out) is an integer used as an error parameter. ier = 33 indicates failure to converge within itmax iterations for at least one of the (n) new roots.

\end{itemize}

683

\end{itemize}

This subroutine always returns the last approximation for root j in x(j). if the convergence criterion is satisfied, then infer(j) is less than or equal to itmax. if the convergence criterion is not satisified, then infer(j) is set to either itmax+1 or itmax+k, with k greater than 1. infer(j) = itmax+1 indicates that muller did not obtain convergence in the allowed number of iterations. in this case, the user may wish to set itmax to a larger value. infer(j) = itmax+k means that convergence was obtained (on iteration k) for the deflated function fp(z) = f(z)/((z-z(1)...(z-z(j-1))) but failed for f(z). in this case, better initial guesses might help or, it might be necessary to relax the convergence criterion.

684

This subroutine always returns the last approximation for root j in x(j). if the convergence criterion is satisfied, then infer(j) is less than or equal to itmax. if the convergence criterion is not satisified, then infer(j) is set to either itmax+1 or itmax+k, with k greater than 1. infer(j) = itmax+1 indicates that muller did not obtain convergence in the allowed number of iterations. in this case, the user may wish to set itmax to a larger value. infer(j) = itmax+k means that convergence was obtained (on iteration k) for the deflated function fp(z) = f(z)/((z-z(1)...(z-z(j-1))) but failed for f(z). in this case, better initial guesses might help or, it might be necessary to relax the convergence criterion.

685

\subsection*{Example}

686

\subsection*{Example}

Example to find the ten roots of $x^{10}-1$

687

Example to find the ten roots of $x^{10}-1$

\begin{verbatim}

688

\begin{verbatim}

program muller

689

program muller

use fvn

690

use fvn

implicit none

691

implicit none

692

integer :: i,info

693

integer :: i,info

complex(8),dimension(10) :: roots

694

complex(8),dimension(10) :: roots

integer,dimension(10) :: infer

695

integer,dimension(10) :: infer

complex(8), external :: f

696

complex(8), external :: f

697

call fvn_z_muller(f,1.d-12,1.d-10,0,0,10,roots,200,infer,info)

698

call fvn_z_muller(f,1.d-12,1.d-10,0,0,10,roots,200,infer,info)

699

write(*,*) "Error code :",info

700

write(*,*) "Error code :",info

do i=1,10

701

do i=1,10

write(*,*) roots(i),infer(i)

702

write(*,*) roots(i),infer(i)

enddo

703

enddo

end program

704

end program

705

function f(x)

706

function f(x)

complex(8) :: x,f

707

complex(8) :: x,f

f=x**10-1

708

f=x**10-1

end function

709

end function

710

\end{verbatim}

711

\end{verbatim}

712

713

714

\section{Numerical integration}

715

\section{Numerical integration}

Using an integrated slightly modified version of quadpack \url{http://www.netlib.org/quadpack}, fvn provide adaptative numerical integration (Gauss Kronrod) of real functions of 1 and 2 variables. fvn also provide a function to calculate Gauss-Legendre abscissas and weight, and a simple non adaptative integration subroutine. All routines exists only in fvn for double precision real.

716

Using an integrated slightly modified version of quadpack \url{http://www.netlib.org/quadpack}, fvn provide adaptative numerical integration (Gauss Kronrod) of real functions of 1 and 2 variables. fvn also provide a function to calculate Gauss-Legendre abscissas and weight, and a simple non adaptative integration subroutine. All routines exists only in fvn for double precision real.

717

\subsection{Gauss Legendre Abscissas and Weigth}

718

\subsection{Gauss Legendre Abscissas and Weigth}

This subroutine was inspired by Numerical Recipes routine gauleg.

719

This subroutine was inspired by Numerical Recipes routine gauleg.

\begin{verbatim}

720

\begin{verbatim}

call fvn_gauss_legendre(n,qx,qw)

721

call fvn_gauss_legendre(n,qx,qw)

\end{verbatim}

722

\end{verbatim}

\begin{itemize}

723

\begin{itemize}

\item n (in) is an integer equal to the number of Gauss Legendre points

724

\item n (in) is an integer equal to the number of Gauss Legendre points

\item qx (out) is a real(8) vector of length n containing the abscissas.

725

\item qx (out) is a real(8) vector of length n containing the abscissas.

\item qw (out) is a real(8) vector of length n containing the weigths.

726

\item qw (out) is a real(8) vector of length n containing the weigths.

\end{itemize}

727

\end{itemize}

This subroutine computes n Gauss-Legendre abscissas and weigths

728

This subroutine computes n Gauss-Legendre abscissas and weigths

729

\subsection{Gauss Legendre Numerical Integration}

730

\subsection{Gauss Legendre Numerical Integration}

\begin{verbatim}

731

\begin{verbatim}

call fvn_gl_integ(f,a,b,n,res)

732

call fvn_gl_integ(f,a,b,n,res)

\end{verbatim}

733

\end{verbatim}

\begin{itemize}

734

\begin{itemize}

\item f (in) is a real(8) function to integrate

735

\item f (in) is a real(8) function to integrate

\item a (in) and b (in) are real(8) respectively lower and higher bound of integration

736

\item a (in) and b (in) are real(8) respectively lower and higher bound of integration

\item n (in) is an integer equal to the number of Gauss Legendre points to use

737

\item n (in) is an integer equal to the number of Gauss Legendre points to use

\item res (out) is a real(8) containing the result

738

\item res (out) is a real(8) containing the result

\end{itemize}

739

\end{itemize}

This function is a simple Gauss Legendre integration subroutine, which evaluate the integral of function f as in equation \ref{intsple} using n Gauss-Legendre pairs.

740

This function is a simple Gauss Legendre integration subroutine, which evaluate the integral of function f as in equation \ref{intsple} using n Gauss-Legendre pairs.

741

\subsection{Gauss Kronrod Adaptative Integration}

742

\subsection{Gauss Kronrod Adaptative Integration}

This kind of numerical integration is an iterative procedure which try to achieve a given precision.

743

This kind of numerical integration is an iterative procedure which try to achieve a given precision.

\subsubsection{Numerical integration of a one variable function}

744

\subsubsection{Numerical integration of a one variable function}

\begin{verbatim}

745

\begin{verbatim}

call fvn_integ_1_gk(f,a,b,epsabs,epsrel,key,res,abserr,ier,limit)

746

call fvn_integ_1_gk(f,a,b,epsabs,epsrel,key,res,abserr,ier,limit)

\end{verbatim}

747

\end{verbatim}

This routine evaluate the integral of function f as in equation \ref{intsple}

748

This routine evaluate the integral of function f as in equation \ref{intsple}

\begin{itemize}

749

\begin{itemize}

\item f (in) is an external real(8) function of one variable

750

\item f (in) is an external real(8) function of one variable

\item a (in) and b (in) are real(8) respectively lower an higher bound of integration

751

\item a (in) and b (in) are real(8) respectively lower an higher bound of integration

\item epsabs (in) and epsrel (in) are real(8) respectively desired absolute and relative error

752

\item epsabs (in) and epsrel (in) are real(8) respectively desired absolute and relative error

\item key (in) is an integer between 1 and 6 correspondind to the Gauss-Kronrod rule to use :

753

\item key (in) is an integer between 1 and 6 correspondind to the Gauss-Kronrod rule to use :

\begin{itemize}

754

\begin{itemize}

\item 1 : 7 - 15 points

755

\item 1 : 7 - 15 points

\item 2 : 10 - 21 points

756

\item 2 : 10 - 21 points

\item 3 : 15 - 31 points

757

\item 3 : 15 - 31 points

\item 4 : 20 - 41 points

758

\item 4 : 20 - 41 points

\item 5 : 25 - 51 points

759

\item 5 : 25 - 51 points

\item 6 : 30 - 61 points

760

\item 6 : 30 - 61 points

\end{itemize}

761

\end{itemize}

\item res (out) is a real(8) containing the estimation of the integration.

762

\item res (out) is a real(8) containing the estimation of the integration.

\item abserr (out) is a real(8) equal to the estimated absolute error

763

\item abserr (out) is a real(8) equal to the estimated absolute error

\item ier (out) is an integer used as an error flag

764

\item ier (out) is an integer used as an error flag

\begin{itemize}

765

\begin{itemize}

\item 0 : no error

766

\item 0 : no error

\item 1 : maximum number of subdivisions allowed has been achieved. one can allow more subdivisions by increasing the value of limit (and taking the according dimension adjustments into account). however, if this yield no improvement it is advised to analyze the integrand in order to determine the integration difficulaties. If the position of a local difficulty can be determined (i.e.singularity, discontinuity within the interval) one will probably gain from splitting up the interval at this point and calling the integrator on the subranges. If possible, an appropriate special-purpose integrator should be used which is designed for handling the type of difficulty involved.

767

\item 1 : maximum number of subdivisions allowed has been achieved. one can allow more subdivisions by increasing the value of limit (and taking the according dimension adjustments into account). however, if this yield no improvement it is advised to analyze the integrand in order to determine the integration difficulaties. If the position of a local difficulty can be determined (i.e.singularity, discontinuity within the interval) one will probably gain from splitting up the interval at this point and calling the integrator on the subranges. If possible, an appropriate special-purpose integrator should be used which is designed for handling the type of difficulty involved.

\item 2 : the occurrence of roundoff error is detected, which prevents the requested tolerance from being achieved.

768

\item 2 : the occurrence of roundoff error is detected, which prevents the requested tolerance from being achieved.

\item 3 : extremely bad integrand behaviour occurs at some points of the integration interval.

769

\item 3 : extremely bad integrand behaviour occurs at some points of the integration interval.

\item 6 : the input is invalid, because (epsabs.le.0 and epsrel.lt.max(50*rel.mach.acc.,0.5d-28)) or limit.lt.1 or lenw.lt.limit*4. result, abserr, neval, last are set to zero. Except when lenw is invalid, iwork(1), work(limit*2+1) and work(limit*3+1) are set to zero, work(1) is set to a and work(limit+1) to b.

770

\item 6 : the input is invalid, because (epsabs.le.0 and epsrel.lt.max(50*rel.mach.acc.,0.5d-28)) or limit.lt.1 or lenw.lt.limit*4. result, abserr, neval, last are set to zero. Except when lenw is invalid, iwork(1), work(limit*2+1) and work(limit*3+1) are set to zero, work(1) is set to a and work(limit+1) to b.

\end{itemize}

771

\end{itemize}

\item limit (in) is an integer equal to maximum number of subintervals in the partition of the given integration interval (a,b). A value of 500 will usually give good results.

772

\item limit (in) is an optional integer equal to maximum number of subintervals in the partition of the given integration interval (a,b). If the parameter is not present a default value of 500 will be used.

\end{itemize}

773

\end{itemize}

774

\begin{equation}

775

\begin{equation}

\int_a^b f(x)~dx

776

\int_a^b f(x)~dx

\label{intsple}

777

\label{intsple}

\end{equation}

778

\end{equation}

779

780

781

782

\subsubsection{Numerical integration of a two variable function}

783

\subsubsection{Numerical integration of a two variable function}

\begin{verbatim}

784

\begin{verbatim}

call fvn_integ_2_gk(f,a,b,g,h,epsabs,epsrel,key,res,abserr,ier,limit)

785

call fvn_integ_2_gk(f,a,b,g,h,epsabs,epsrel,key,res,abserr,ier,limit)

\end{verbatim}

786

\end{verbatim}

This function evaluate the integral of a function f(x,y) as defined in equation \ref{intdble}. The parameters of same name as in the previous paragraph have exactly the same function and behaviour thus only what differs is decribed here

787

This function evaluate the integral of a function f(x,y) as defined in equation \ref{intdble}. The parameters of same name as in the previous paragraph have exactly the same function and behaviour thus only what differs is decribed here

\begin{itemize}

788

\begin{itemize}

\item a (in) and b (in) are real(8) corresponding respectively to lower and higher bound of integration for the x variable.

789

\item a (in) and b (in) are real(8) corresponding respectively to lower and higher bound of integration for the x variable.

\item g(x) (in) and h(x) (in) are external functions describing the lower and higher bound of integration for the y variable as a function of x.

790

\item g(x) (in) and h(x) (in) are external functions describing the lower and higher bound of integration for the y variable as a function of x.

\end{itemize}

791

\end{itemize}

792

\begin{equation}

793

\begin{equation}

\int_a^b \int_{g(x)}^{h(x)} f(x,y)~dy~dx

794

\int_a^b \int_{g(x)}^{h(x)} f(x,y)~dy~dx

\label{intdble}

795

\label{intdble}

\end{equation}

796

\end{equation}

797

\subsubsection*{Example}

798

\subsubsection*{Example}

\begin{verbatim}

799

\begin{verbatim}

program integ

800

program integ

use fvn

801

use fvn

implicit none

802

implicit none

803

real(8), external :: f1,f2,g,h

804

real(8), external :: f1,f2,g,h

real(8) :: a,b,epsabs,epsrel,abserr,res

805

real(8) :: a,b,epsabs,epsrel,abserr,res

integer :: key,ier

806

integer :: key,ier

807

a=0.

808

a=0.

b=1.

809

b=1.

epsabs=1d-8

810

epsabs=1d-8

epsrel=1d-8

811

epsrel=1d-8

key=2

812

key=2

call fvn_d_integ_1_gk(f1,a,b,epsabs,epsrel,key,res,abserr,ier,500)

813

call fvn_d_integ_1_gk(f1,a,b,epsabs,epsrel,key,res,abserr,ier,500)

write(*,*) "Integration of x*x between 0 and 1 : "

814

write(*,*) "Integration of x*x between 0 and 1 : "

write(*,*) res

815

write(*,*) res

816

call fvn_d_integ_2_gk(f2,a,b,g,h,epsabs,epsrel,key,res,abserr,ier,500)

817

call fvn_d_integ_2_gk(f2,a,b,g,h,epsabs,epsrel,key,res,abserr,ier,500)

write(*,*) "Integration of x*y between 0 and 1 on both x and y : "

818

write(*,*) "Integration of x*y between 0 and 1 on both x and y : "

write(*,*) res

819

write(*,*) res

820

821

end program

822

end program

823

function f1(x)

824

function f1(x)

implicit none

825

implicit none

real(8) :: x,f1

826

real(8) :: x,f1

f1=x*x

827

f1=x*x

end function

828

end function

829

function f2(x,y)

830

function f2(x,y)

implicit none

831

implicit none

real(8) :: x,y,f2

832

real(8) :: x,y,f2

f2=x*y

833

f2=x*y

end function

834

end function

835

function g(x)

836

function g(x)

implicit none

837

implicit none

real(8) :: x,g

838

real(8) :: x,g

g=0.

839

g=0.

end function

840

end function

841

function h(x)

842

function h(x)

implicit none

843

implicit none

real(8) :: x,h

844

real(8) :: x,h

h=1.

845

h=1.

end function

846

end function

\end{verbatim}

847

\end{verbatim}

848

849

\section{Special functions}

850

\section{Special functions}

Specials functions are available in fvn by using an implementation of fnlib \url{http://www.netlib.org/fn}. This can be used separatly from the rest of fvn by using the module \verb'fvn_fnlib' and linking the library \verb'libfvn_fnlib.a' . The module provides a generic interfaces to all the routines. Specific names of the routines are given in the description. The double complex versions of the routines are not present in the web version of fnlib, so these have been added, but not intensely tested.

851

Specials functions are available in fvn by using an implementation of fnlib \url{http://www.netlib.org/fn}. This can be used separatly from the rest of fvn by using the module \verb'fvn_fnlib' and linking the library \verb'libfvn_fnlib.a' . The module provides a generic interfaces to all the routines. Specific names of the routines are given in the description. The double complex versions of the routines are not present in the web version of fnlib, so these have been added, but not intensely tested.

852

\paragraph{Important Note}

853

\paragraph{Important Note}

Due to the addition of fnlib to fvn, some functions that were in fvn and are redondant will be removed from fvn, so update your code now and replace them with the fnlib version. These are listed here after :

854

Due to the addition of fnlib to fvn, some functions that were in fvn and are redondant will be removed from fvn, so update your code now and replace them with the fnlib version. These are listed here after :

\begin{itemize}

855

\begin{itemize}

\item \verb'fvn_z_acos' replaced by \verb'acos'

856

\item \verb'fvn_z_acos' replaced by \verb'acos'

\item \verb'fvn_z_asin' replaced by \verb'asin'

857

\item \verb'fvn_z_asin' replaced by \verb'asin'

\item \verb'fvn_d_asinh' replaced by \verb'asinh'

858

\item \verb'fvn_d_asinh' replaced by \verb'asinh'

\item \verb'fvn_d_acosh' replaced by \verb'acosh'

859

\item \verb'fvn_d_acosh' replaced by \verb'acosh'

\item \verb'fvn_s_csevl' replaced by \verb'csevl'

860

\item \verb'fvn_s_csevl' replaced by \verb'csevl'

\item \verb'fvn_d_csevl' replaced by \verb'csevl'

861

\item \verb'fvn_d_csevl' replaced by \verb'csevl'

\item \verb'fvn_d_factorial' replaced by \verb'fac'

862

\item \verb'fvn_d_factorial' replaced by \verb'fac'

\item \verb'fvn_d_lngamma' replaced by \verb'alngam'

863

\item \verb'fvn_d_lngamma' replaced by \verb'alngam'

\end{itemize}

864

\end{itemize}

865

866

\subsection{Elementary functions}

867

\subsection{Elementary functions}

\subsubsection{carg}

868

\subsubsection{carg}

\begin{verbatim}

869

\begin{verbatim}

carg(z)

870

carg(z)

\end{verbatim}

871

\end{verbatim}

\begin{itemize}

872

\begin{itemize}

\item z (in) is a complex

873

\item z (in) is a complex

\end{itemize}

874

\end{itemize}

This function evaluates the argument of the complex z. That is $\theta$ for $z=\rho e^{i\theta}$.

875

This function evaluates the argument of the complex z. That is $\theta$ for $z=\rho e^{i\theta}$.

876

Specific interfaces : \verb'carg,zarg'

877

Specific interfaces : \verb'carg,zarg'

878

879

\subsubsection{cbrt}

880

\subsubsection{cbrt}

\begin{verbatim}

881

\begin{verbatim}

cbrt(x)

882

cbrt(x)

\end{verbatim}

883

\end{verbatim}

\begin{itemize}

884

\begin{itemize}

\item x is a real or complex

885

\item x is a real or complex

\end{itemize}

886

\end{itemize}

This function evaluates the cubic root of the argument x.

887

This function evaluates the cubic root of the argument x.

888

Specific interfaces : \verb'cbrt,dcbrt,ccbrt,zcbrt'

889

Specific interfaces : \verb'cbrt,dcbrt,ccbrt,zcbrt'

890

891

\subsubsection{exprl}

892

\subsubsection{exprl}

\begin{verbatim}

893

\begin{verbatim}

exprl(x)

894

exprl(x)

\end{verbatim}

895

\end{verbatim}

\begin{itemize}

896

\begin{itemize}

\item x is a real or complex

897

\item x is a real or complex

\end{itemize}

898

\end{itemize}

This function evaluates ${e^x-1}\over x$.

899

This function evaluates ${e^x-1}\over x$.

900

Specific interfaces : \verb'exprel,dexprl,cexprl,zexprl'

901

Specific interfaces : \verb'exprel,dexprl,cexprl,zexprl'

902

\subsubsection{log10}

903

\subsubsection{log10}

\begin{verbatim}

904

\begin{verbatim}

log10(x)

905

log10(x)

\end{verbatim}

906

\end{verbatim}

\begin{itemize}

907

\begin{itemize}

\item x is a real or complex

908

\item x is a real or complex

\end{itemize}

909

\end{itemize}

This function is an extension of the intrinsic function log10 to complex arguments.

910

This function is an extension of the intrinsic function log10 to complex arguments.

911

Specific interfaces : \verb'clog10,zlog10'

912

Specific interfaces : \verb'clog10,zlog10'

913

914

\subsubsection{alnrel}

915

\subsubsection{alnrel}

\begin{verbatim}

916

\begin{verbatim}

alnrel(x)

917

alnrel(x)

\end{verbatim}

918

\end{verbatim}

\begin{itemize}

919

\begin{itemize}

\item x is a real or complex

920

\item x is a real or complex

\end{itemize}

921

\end{itemize}

This function evaluates $ln(1+x)$.

922

This function evaluates $ln(1+x)$.

923

Specific interfaces : \verb'alnrel,dlnrel,clnrel,zlnrel'

924

Specific interfaces : \verb'alnrel,dlnrel,clnrel,zlnrel'

925

926

\subsection{Trigonometry}

927

\subsection{Trigonometry}

\subsubsection{tan}

928

\subsubsection{tan}

\begin{verbatim}

929

\begin{verbatim}

tan(x)

930

tan(x)

\end{verbatim}

931

\end{verbatim}

\begin{itemize}

932

\begin{itemize}

\item x is a real or complex

933

\item x is a real or complex

\end{itemize}

934

\end{itemize}

This function evaluates the tangent of the argument. It is an extension of the intrinsic function tan to complex arguments.

935

This function evaluates the tangent of the argument. It is an extension of the intrinsic function tan to complex arguments.

936

Specific interfaces : \verb'ctan,ztan'

937

Specific interfaces : \verb'ctan,ztan'

938

939

\subsubsection{cot}

940

\subsubsection{cot}

\begin{verbatim}

941

\begin{verbatim}

cot(x)

942

cot(x)

\end{verbatim}

943

\end{verbatim}

\begin{itemize}

944

\begin{itemize}

\item x is a real or complex

945

\item x is a real or complex

\end{itemize}

946

\end{itemize}

This function evaluate the cotangent of the argument.

947

This function evaluate the cotangent of the argument.

948

Specific interfaces : \verb'cot,dcot,ccot,zcot'

949

Specific interfaces : \verb'cot,dcot,ccot,zcot'

950

\subsubsection{sindg}

951

\subsubsection{sindg}

\begin{verbatim}

952

\begin{verbatim}

sindg(x)

953

sindg(x)

\end{verbatim}

954

\end{verbatim}

\begin{itemize}

955

\begin{itemize}

\item x is a real

956

\item x is a real

\end{itemize}

957

\end{itemize}

This function evaluate the sinus of the argument expressed in degrees.

958

This function evaluate the sinus of the argument expressed in degrees.

959

Specific interfaces : \verb'sindg,dsindg'

960

Specific interfaces : \verb'sindg,dsindg'

961

962

\subsubsection{cosdg}

963

\subsubsection{cosdg}

\begin{verbatim}

964

\begin{verbatim}

cosdg(x)

965

cosdg(x)

\end{verbatim}

966

\end{verbatim}

\begin{itemize}

967

\begin{itemize}

\item x is a real

968

\item x is a real

\end{itemize}

969

\end{itemize}

This function evaluate the cosinus of the argument expressed in degrees.

970

This function evaluate the cosinus of the argument expressed in degrees.

971

Specific interfaces : \verb'cosdg,dcosdg'

972

Specific interfaces : \verb'cosdg,dcosdg'

973

974

\subsubsection{asin}

975

\subsubsection{asin}

\begin{verbatim}

976

\begin{verbatim}

asin(x)

977

asin(x)

\end{verbatim}

978

\end{verbatim}

\begin{itemize}

979

\begin{itemize}

\item x is a real or complex

980

\item x is a real or complex

\end{itemize}

981

\end{itemize}

This function evaluates the arc sine of the argument. It is an extension of the intrinsic function asin to complex arguments.

982

This function evaluates the arc sine of the argument. It is an extension of the intrinsic function asin to complex arguments.

983

Specific interfaces : \verb'casin,zasin'

984

Specific interfaces : \verb'casin,zasin'

985

\subsubsection{acos}

986

\subsubsection{acos}

\begin{verbatim}

987

\begin{verbatim}

acos(x)

988

acos(x)

\end{verbatim}

989

\end{verbatim}

\begin{itemize}

990

\begin{itemize}

\item x is a real or complex

991

\item x is a real or complex

\end{itemize}

992

\end{itemize}

This function evaluates the arc cosine of the argument. It is an extension of the intrinsic function acos to complex arguments.

993

This function evaluates the arc cosine of the argument. It is an extension of the intrinsic function acos to complex arguments.

994

Specific interfaces : \verb'cacos,zacos'

995

Specific interfaces : \verb'cacos,zacos'

996

997

\subsubsection{atan}

998

\subsubsection{atan}

\begin{verbatim}

999

\begin{verbatim}

atan(x)

1000

atan(x)

\end{verbatim}

1001

\end{verbatim}

\begin{itemize}

1002

\begin{itemize}

\item x is a real or complex

1003

\item x is a real or complex

\end{itemize}

1004

\end{itemize}

This function evaluates the arc tangent of the argument. It is an extension of the intrinsic function atan to complex arguments.

1005

This function evaluates the arc tangent of the argument. It is an extension of the intrinsic function atan to complex arguments.

1006

Specific interfaces : \verb'catan,zatan'

1007

Specific interfaces : \verb'catan,zatan'

1008

\subsubsection{atan2}

1009

\subsubsection{atan2}

\begin{verbatim}

1010

\begin{verbatim}

atan2(x,y)

1011

atan2(x,y)

\end{verbatim}

1012

\end{verbatim}

\begin{itemize}

1013

\begin{itemize}

\item x,y are real or complex

1014

\item x,y are real or complex

\end{itemize}

1015

\end{itemize}

This function evaluates the arc tangent of $x \over y$. It is an extension of the intrinsic function atan2 to complex arguments.

1016

This function evaluates the arc tangent of $x \over y$. It is an extension of the intrinsic function atan2 to complex arguments.

1017

Specific interfaces : \verb'catan2,zatan2'

1018

Specific interfaces : \verb'catan2,zatan2'

1019

\subsubsection{sinh}

1020

\subsubsection{sinh}

\begin{verbatim}

1021

\begin{verbatim}

sinh(x)

1022

sinh(x)

\end{verbatim}

1023

\end{verbatim}

\begin{itemize}

1024

\begin{itemize}

\item x is a real or complex

1025

\item x is a real or complex

\end{itemize}

1026

\end{itemize}

This function evaluates the hyperbolic sine of the argument. It is an extension of the intrinsic function sinh to complex arguments.

1027

This function evaluates the hyperbolic sine of the argument. It is an extension of the intrinsic function sinh to complex arguments.

1028

Specific interfaces : \verb'csinh,zsinh'

1029

Specific interfaces : \verb'csinh,zsinh'

1030

1031

\subsubsection{cosh}

1032

\subsubsection{cosh}

\begin{verbatim}

1033

\begin{verbatim}

cosh(x)

1034

cosh(x)

\end{verbatim}

1035

\end{verbatim}

\begin{itemize}

1036

\begin{itemize}

\item x is a real or complex

1037

\item x is a real or complex

\end{itemize}

1038

\end{itemize}

This function evaluates the hyperbolic cosine of the argument. It is an extension of the intrinsic function cosh to complex arguments.

1039

This function evaluates the hyperbolic cosine of the argument. It is an extension of the intrinsic function cosh to complex arguments.

1040

Specific interfaces : \verb'ccosh,zcosh'

1041

Specific interfaces : \verb'ccosh,zcosh'

1042

\subsubsection{tanh}

1043

\subsubsection{tanh}

\begin{verbatim}

1044

\begin{verbatim}

tanh(x)

1045

tanh(x)

\end{verbatim}

1046

\end{verbatim}

This function evaluates the hyperbolic tangent of the argument. It is an extension of the intrinsic function tanh to complex arguments.

1047

This function evaluates the hyperbolic tangent of the argument. It is an extension of the intrinsic function tanh to complex arguments.

1048

Specific interfaces : \verb'ctanh,ztanh'

1049

Specific interfaces : \verb'ctanh,ztanh'

1050

\subsubsection{asinh}

1051

\subsubsection{asinh}

\begin{verbatim}

1052

\begin{verbatim}

asinh(x)

1053

asinh(x)

\end{verbatim}

1054

\end{verbatim}

\begin{itemize}

1055

\begin{itemize}

\item x is a real or complex

1056

\item x is a real or complex

\end{itemize}

1057

\end{itemize}

This function evaluates the arc hyperbolic sine of the argument.

1058

This function evaluates the arc hyperbolic sine of the argument.

1059

Specific interfaces : \verb'asinh,dasinh,casinh,zasinh'

1060

Specific interfaces : \verb'asinh,dasinh,casinh,zasinh'

1061

\subsubsection{acosh}

1062

\subsubsection{acosh}

\begin{verbatim}

1063

\begin{verbatim}

acosh(x)

1064

acosh(x)

\end{verbatim}

1065

\end{verbatim}

\begin{itemize}

1066

\begin{itemize}

\item x is a real or complex

1067

\item x is a real or complex

\end{itemize}

1068

\end{itemize}

This function evaluates the arc hyperbolic cosine of the argument.

1069

This function evaluates the arc hyperbolic cosine of the argument.

1070

Specific interfaces : \verb'acosh,dacosh,cacosh,zacosh'

1071

Specific interfaces : \verb'acosh,dacosh,cacosh,zacosh'

1072

\subsubsection{atanh}

1073

\subsubsection{atanh}

\begin{verbatim}

1074

\begin{verbatim}

atanh(x)

1075

atanh(x)

\end{verbatim}

1076

\end{verbatim}

\begin{itemize}

1077

\begin{itemize}

\item x is a real or complex

1078

\item x is a real or complex

\end{itemize}

1079

\end{itemize}

This function evaluates the arc hyperbolic tangent of the argument.

1080

This function evaluates the arc hyperbolic tangent of the argument.

1081

Specific interfaces : \verb'atanh,datanh,catanh,zatanh'

1082

Specific interfaces : \verb'atanh,datanh,catanh,zatanh'

1083

\subsection{Exponential Integral and related}

1084

\subsection{Exponential Integral and related}

\subsubsection{ei}

1085

\subsubsection{ei}

\begin{verbatim}

1086

\begin{verbatim}

ei(x)

1087

ei(x)

\end{verbatim}

1088

\end{verbatim}

\begin{itemize}

1089

\begin{itemize}

\item x is a real

1090

\item x is a real

\end{itemize}

1091

\end{itemize}

This function evaluates the exponential integral for argument greater then 0 and the Cauchy principal value for argument less than 0. It is define by equation \ref{ei} for $x \neq 0$.

1092

This function evaluates the exponential integral for argument greater then 0 and the Cauchy principal value for argument less than 0. It is define by equation \ref{ei} for $x \neq 0$.

\begin{equation}

1093

\begin{equation}

\label{ei}

1094

\label{ei}

ei(x)= - \int _{-x} ^\infty {e^{-t}\over t}dt

1095

ei(x)= - \int _{-x} ^\infty {e^{-t}\over t}dt

\end{equation}

1096

\end{equation}

1097

Specific interfaces : \verb'ei,dei'

1098

Specific interfaces : \verb'ei,dei'

1099

1100

\subsubsection{e1}

1101

\subsubsection{e1}

\begin{verbatim}

1102

\begin{verbatim}

e1(x)

1103

e1(x)

\end{verbatim}

1104

\end{verbatim}

\begin{itemize}

1105

\begin{itemize}

\item x is a real

1106

\item x is a real

\end{itemize}

1107

\end{itemize}

This function evaluates the exponential integral for argument greater than 0 and the Cauchy principal value for argument less than 0. It is define by equation \ref{e1} for $x \neq 0$.

1108

This function evaluates the exponential integral for argument greater than 0 and the Cauchy principal value for argument less than 0. It is define by equation \ref{e1} for $x \neq 0$.

\begin{equation}

1109

\begin{equation}

\label{e1}

1110

\label{e1}

e1(x)= \int _{x} ^\infty {e^{-t}\over t}dt

1111

e1(x)= \int _{x} ^\infty {e^{-t}\over t}dt

\end{equation}

1112

\end{equation}

1113

Specific interfaces : \verb'e1,de1'

1114

Specific interfaces : \verb'e1,de1'

1115

\subsubsection{ali}

1116

\subsubsection{ali}

\begin{verbatim}

1117

\begin{verbatim}

ali(x)

1118

ali(x)

\end{verbatim}

1119

\end{verbatim}

\begin{itemize}

1120

\begin{itemize}

\item x is a real

1121

\item x is a real

\end{itemize}

1122

\end{itemize}

This function evaluates the logarithm integral. it is define by equation \ref{ali} for $x > 0$ and $x \neq 1$.

1123

This function evaluates the logarithm integral. it is define by equation \ref{ali} for $x > 0$ and $x \neq 1$.

\begin{equation}

1124

\begin{equation}

\label{ali}

1125

\label{ali}

ali(x)= - \int _0 ^x {dt \over ln(x)}

1126

ali(x)= - \int _0 ^x {dt \over ln(x)}

\end{equation}

1127

\end{equation}

1128

Specific interfaces : \verb'ali,dli'

1129

Specific interfaces : \verb'ali,dli'

1130

\subsubsection{si}

1131

\subsubsection{si}

\begin{verbatim}

1132

\begin{verbatim}

si(x)

1133

si(x)

\end{verbatim}

1134

\end{verbatim}

\begin{itemize}

1135

\begin{itemize}

\item x is a real

1136

\item x is a real

\end{itemize}

1137

\end{itemize}

This function evaluates the sine integral defined by equation \ref{si}.

1138

This function evaluates the sine integral defined by equation \ref{si}.

\begin{equation}

1139

\begin{equation}

\label{si}

1140

\label{si}

si(x)= \int _0 ^x {sin(t) \over t }dt

1141

si(x)= \int _0 ^x {sin(t) \over t }dt

\end{equation}

1142

\end{equation}

1143

Specific interfaces : \verb'si,dsi'

1144

Specific interfaces : \verb'si,dsi'

1145

1146

\subsubsection{ci}

1147

\subsubsection{ci}

\begin{verbatim}

1148

\begin{verbatim}

ci(x)

1149

ci(x)

\end{verbatim}

1150

\end{verbatim}

\begin{itemize}

1151

\begin{itemize}

\item x is a real

1152

\item x is a real

\end{itemize}

1153

\end{itemize}

This function evaluates the cosine integral defined by equation \ref{ci} where $\gamma \approx 0.57721566$ represent Euler's constant.

1154

This function evaluates the cosine integral defined by equation \ref{ci} where $\gamma \approx 0.57721566$ represent Euler's constant.

\begin{equation}

1155

\begin{equation}

\label{ci}

1156

\label{ci}

ci(x)= \gamma + ln(x) + \int _0 ^x {{1-cos(t)} \over t} dt

1157

ci(x)= \gamma + ln(x) + \int _0 ^x {{1-cos(t)} \over t} dt

\end{equation}

1158

\end{equation}

1159

Specific interfaces : \verb'ci,dci'

1160

Specific interfaces : \verb'ci,dci'

1161

\subsubsection{cin}

1162

\subsubsection{cin}

\begin{verbatim}

1163

\begin{verbatim}

cin(x)

1164

cin(x)

\end{verbatim}

1165

\end{verbatim}

\begin{itemize}

1166

\begin{itemize}

\item x is a real

1167

\item x is a real

\end{itemize}

1168

\end{itemize}

This function evaluates the cosine integral alternate definition given by equation \ref{cin}.

1169

This function evaluates the cosine integral alternate definition given by equation \ref{cin}.

\begin{equation}

1170

\begin{equation}

\label{cin}

1171

\label{cin}

cin(x)= \int _0 ^x {{1-cos(t)} \over t} dt

1172

cin(x)= \int _0 ^x {{1-cos(t)} \over t} dt

\end{equation}

1173

\end{equation}

1174

Specific interface : \verb'cin,dcin'

1175

Specific interface : \verb'cin,dcin'

1176

\subsubsection{shi}

1177

\subsubsection{shi}

\begin{equation}

1178

\begin{equation}

shi(x)

1179

shi(x)

\end{equation}

1180

\end{equation}

\begin{itemize}

1181

\begin{itemize}

\item x is a real

1182

\item x is a real

\end{itemize}

1183

\end{itemize}

This function evaluates the hyperbolic sine integral defined by equation \ref{shi}.

1184

This function evaluates the hyperbolic sine integral defined by equation \ref{shi}.

\begin{equation}

1185

\begin{equation}

\label{shi}

1186

\label{shi}

shi(x) = \int _0 ^x {sinh(t) \over t}dt

1187

shi(x) = \int _0 ^x {sinh(t) \over t}dt

\end{equation}

1188

\end{equation}

1189

Specific interfaces : \verb'shi,dshi'

1190

Specific interfaces : \verb'shi,dshi'

1191

\subsubsection{chi}

1192

\subsubsection{chi}

\begin{verbatim}

1193

\begin{verbatim}

chi(x)

1194

chi(x)

\end{verbatim}

1195

\end{verbatim}

\begin{itemize}

1196

\begin{itemize}

\item x is a real

1197

\item x is a real

\end{itemize}

1198

\end{itemize}

This function evaluates the hyperbolic cosine integral defined by equation \ref{chi} where $\gamma \approx 0.57721566$ represent Euler's constant.

1199

This function evaluates the hyperbolic cosine integral defined by equation \ref{chi} where $\gamma \approx 0.57721566$ represent Euler's constant.

\begin{equation}

1200

\begin{equation}

\label{chi}

1201

\label{chi}

chi(x)= \gamma + ln(x) + \int _0 ^x {{cosh(t) -1} \over t}dt

1202

chi(x)= \gamma + ln(x) + \int _0 ^x {{cosh(t) -1} \over t}dt

\end{equation}

1203

\end{equation}

1204

Specific interfaces : chi,dchi

1205

Specific interfaces : chi,dchi

1206

1207

\subsubsection{cinh}

1208

\subsubsection{cinh}

\begin{verbatim}

1209

\begin{verbatim}

cinh(x)

1210

cinh(x)

\end{verbatim}

1211

\end{verbatim}

\begin{itemize}

1212

\begin{itemize}

\item x is a real

1213

\item x is a real

\end{itemize}

1214

\end{itemize}

This function evaluates the hyperbolic cosine integral alternate definition given by equation \ref{cinh}.

1215

This function evaluates the hyperbolic cosine integral alternate definition given by equation \ref{cinh}.

\begin{equation}

1216

\begin{equation}

\label{cinh}

1217

\label{cinh}

cinh(x) = \int _0 ^x {{cosh(t) -1} \over t}dt

1218

cinh(x) = \int _0 ^x {{cosh(t) -1} \over t}dt

\end{equation}

1219

\end{equation}

1220

Specific interfaces : cinh,dcinh

1221

Specific interfaces : cinh,dcinh

1222

1223

\subsection{Gamma function and related}

1224

\subsection{Gamma function and related}

\subsubsection{fac}

1225

\subsubsection{fac}

\begin{verbatim}

1226

\begin{verbatim}

fac(n)

1227

fac(n)

dfac(n)

1228

dfac(n)

\end{verbatim}

1229

\end{verbatim}

\begin{itemize}

1230

\begin{itemize}

\item n is an integer

1231

\item n is an integer

\end{itemize}

1232

\end{itemize}

This function return $n!$ as a real(4) or real(8) for dfac. There's no generic interface for this one.

1233

This function return $n!$ as a real(4) or real(8) for dfac. There's no generic interface for this one.

1234

Specific interfaces : \verb'fac,dfac'

1235

Specific interfaces : \verb'fac,dfac'

1236

\subsubsection{binom}

1237

\subsubsection{binom}

\begin{verbatim}

1238

\begin{verbatim}

binom(n,m)

1239

binom(n,m)

dbinom(n,m)

1240

dbinom(n,m)

\end{verbatim}

1241

\end{verbatim}

\begin{itemize}

1242

\begin{itemize}

\item n,m are integers

1243

\item n,m are integers

\end{itemize}

1244

\end{itemize}

This function return the binomial coefficient defined by equation \ref{binom} with $n \geq m \geq 0$. binom returns a real(4), dbinom a real(8). There's no generic interface for this one.

1245

This function return the binomial coefficient defined by equation \ref{binom} with $n \geq m \geq 0$. binom returns a real(4), dbinom a real(8). There's no generic interface for this one.

\begin{equation}

1246

\begin{equation}

\label{binom}

1247

\label{binom}

binom(n,m) = C_n^m = {{n!} \over {m!(n-m)!}}

1248

binom(n,m) = C_n^m = {{n!} \over {m!(n-m)!}}

\end{equation}

1249

\end{equation}

1250

Specific interfaces : \verb'binom,dbinom'

1251

Specific interfaces : \verb'binom,dbinom'

1252

1253

\subsubsection{gamma}

1254

\subsubsection{gamma}

\begin{verbatim}

1255

\begin{verbatim}

gamma(x)

1256

gamma(x)

\end{verbatim}

1257

\end{verbatim}

\begin{itemize}

1258

\begin{itemize}

\item x is a real or complex

1259

\item x is a real or complex

\end{itemize}

1260

\end{itemize}

This function evaluates $ \Gamma (x) $ defined by equation \ref{gamma}.

1261

This function evaluates $ \Gamma (x) $ defined by equation \ref{gamma}.

\begin{equation}

1262

\begin{equation}

\label{gamma}

1263

\label{gamma}

\Gamma (x) = \int _0 ^{\infty} t^{x-1}e^{-t}dt

1264

\Gamma (x) = \int _0 ^{\infty} t^{x-1}e^{-t}dt

\end{equation}

1265

\end{equation}

Note that $n!=\Gamma (n+1)$.

1266

Note that $n!=\Gamma (n+1)$.

1267

Specific interfaces :\verb'gamma,dgamma,cgamma,zgamm'

1268

Specific interfaces :\verb'gamma,dgamma,cgamma,zgamm'

1269

\subsubsection{gamr}

1270

\subsubsection{gamr}

\begin{verbatim}

1271

\begin{verbatim}

gamr(x)

1272

gamr(x)

\end{verbatim}

1273

\end{verbatim}

\begin{itemize}

1274

\begin{itemize}

\item x is a real or complex

1275

\item x is a real or complex

\end{itemize}

1276

\end{itemize}

This function evaluates the reciprocal gamma function $gamr(x)= {1 \over \Gamma(x)}$

1277

This function evaluates the reciprocal gamma function $gamr(x)= {1 \over \Gamma(x)}$

1278

1279

\subsubsection{alngam}

1280

\subsubsection{alngam}

\begin{verbatim}

1281

\begin{verbatim}

alngam(x)

1282

alngam(x)

\end{verbatim}

1283

\end{verbatim}

\begin{itemize}

1284

\begin{itemize}

\item x is a real or complex

1285

\item x is a real or complex

\end{itemize}

1286

\end{itemize}

This function evaluates $ln(|\Gamma(x)|)$

1287

This function evaluates $ln(|\Gamma(x)|)$

1288

Specific interfaces : \verb'alngam,dlngam,clngam,zlngam'

1289

Specific interfaces : \verb'alngam,dlngam,clngam,zlngam'

1290

1291

\subsubsection{algams}

1292

\subsubsection{algams}

\begin{verbatim}

1293

\begin{verbatim}

call algams(x,algam,sgngam)

1294

call algams(x,algam,sgngam)

\end{verbatim}

1295

\end{verbatim}

\begin{itemize}

1296

\begin{itemize}

\item x (in) is a real

1297

\item x (in) is a real

\item algam (out) is a real

1298

\item algam (out) is a real

\item sgngam (out) is a real

1299

\item sgngam (out) is a real

\end{itemize}

1300

\end{itemize}

This subroutine evaluates the logarithm of the absolute value of gamma and the sign of gamma.

1301

This subroutine evaluates the logarithm of the absolute value of gamma and the sign of gamma.

$algam=ln(|\Gamma(x)|)$ and $sgngam=1.0$ or $-1.0$ according to the sign of $\Gamma(x)$.

1302

$algam=ln(|\Gamma(x)|)$ and $sgngam=1.0$ or $-1.0$ according to the sign of $\Gamma(x)$.

1303

Specific interfaces : \verb'algams,dlgams'

1304

Specific interfaces : \verb'algams,dlgams'

1305

\subsubsection{gami}

1306

\subsubsection{gami}

\begin{verbatim}

1307

\begin{verbatim}

gami(a,x)

1308

gami(a,x)

\end{verbatim}

1309

\end{verbatim}

\begin{itemize}

1310

\begin{itemize}

\item x is a positive real

1311

\item x is a positive real

\item a is a strictly positive real

1312

\item a is a strictly positive real

\end{itemize}

1313

\end{itemize}

This function evaluates the incomplete gamma function defined by equation \ref{gami}.

1314

This function evaluates the incomplete gamma function defined by equation \ref{gami}.

\begin{equation}

1315

\begin{equation}

\label{gami}

1316

\label{gami}

gami(a,x)=\gamma(a,x)=\int _0 ^x t^{a-1} e^{-t}dt

1317

gami(a,x)=\gamma(a,x)=\int _0 ^x t^{a-1} e^{-t}dt

\end{equation}

1318

\end{equation}

1319

Specific interfaces : \verb'gami,dgami'

1320

Specific interfaces : \verb'gami,dgami'

1321

\subsubsection{gamic}

1322

\subsubsection{gamic}

\begin{verbatim}

1323

\begin{verbatim}

gamic(a,x)

1324

gamic(a,x)

\end{verbatim}

1325

\end{verbatim}

\begin{itemize}

1326

\begin{itemize}

\item x is a positive real

1327

\item x is a positive real

\item a is a real

1328

\item a is a real

\end{itemize}

1329

\end{itemize}

This function evaluates the complementary incomplete gamma function defined by equation \ref{gamic}.

1330

This function evaluates the complementary incomplete gamma function defined by equation \ref{gamic}.

\begin{equation}

1331

\begin{equation}

\label{gamic}

1332

\label{gamic}

gamic(a,x)=\Gamma(a,x)=\int _x ^\infty t^{a-1} e^{-t}dt

1333

gamic(a,x)=\Gamma(a,x)=\int _x ^\infty t^{a-1} e^{-t}dt

\end{equation}

1334

\end{equation}

1335

Specific interfaces : \verb'gamic,dgamic'

1336

Specific interfaces : \verb'gamic,dgamic'

1337

\subsubsection{gamit}

1338

\subsubsection{gamit}

\begin{verbatim}

1339

\begin{verbatim}

gamit(a,x)

1340

gamit(a,x)

\end{verbatim}

1341

\end{verbatim}

\begin{itemize}

1342

\begin{itemize}

\item x is a positive real

1343

\item x is a positive real

\item a is a real

1344

\item a is a real

\end{itemize}

1345

\end{itemize}

This function evaluates the Tricomi's incomplete gamma function defined by equation \ref{gamit}.

1346

This function evaluates the Tricomi's incomplete gamma function defined by equation \ref{gamit}.

\begin{equation}

1347

\begin{equation}

\label{gamit}

1348

\label{gamit}

gamit(a,x)=\gamma^* (a,x)= {{x^{-a}\gamma(a,x)}\over \Gamma(a)}

1349

gamit(a,x)=\gamma^* (a,x)= {{x^{-a}\gamma(a,x)}\over \Gamma(a)}

\end{equation}

1350

\end{equation}

1351

Specific interfaces : \verb'gamit,dgamit'

1352

Specific interfaces : \verb'gamit,dgamit'

1353

1354

\subsubsection{psi}

1355

\subsubsection{psi}

\begin{verbatim}

1356

\begin{verbatim}

psi(x)

1357

psi(x)

\end{verbatim}

1358

\end{verbatim}

\begin{itemize}

1359

\begin{itemize}

\item x is a real or complex

1360

\item x is a real or complex

\end{itemize}

1361

\end{itemize}

This function evaluates the psi function which is the logarithm derivative of the gamma function as defined in equation \ref{psi}.

1362

This function evaluates the psi function which is the logarithm derivative of the gamma function as defined in equation \ref{psi}.

\begin{equation}

1363

\begin{equation}

\label{psi}

1364

\label{psi}

psi(x)= \psi(x) = {d\over dx} ln(\Gamma(x))

1365

psi(x)= \psi(x) = {d\over dx} ln(\Gamma(x))

\end{equation}

1366

\end{equation}

x must not be zero or a negative integer.

1367

x must not be zero or a negative integer.

1368

Specific interfaces : \verb'psi,dpsi,cpsi,zpsi'

1369

Specific interfaces : \verb'psi,dpsi,cpsi,zpsi'

1370

1371

\subsubsection{poch}

1372

\subsubsection{poch}

\begin{verbatim}

1373

\begin{verbatim}

poch(a,x)

1374

poch(a,x)

\end{verbatim}

1375

\end{verbatim}

\begin{itemize}

1376

\begin{itemize}

\item x is a real

1377

\item x is a real

\item a is a real

1378

\item a is a real

\end{itemize}

1379

\end{itemize}

This function evaluates a generalization of Pochhammer's symbol.

1380

This function evaluates a generalization of Pochhammer's symbol.

1381

Pochhammer's symbol for n a positive integer is given by equation \ref{poch_int}

1382

Pochhammer's symbol for n a positive integer is given by equation \ref{poch_int}

\begin{equation}

1383

\begin{equation}

\label{poch_int}

1384

\label{poch_int}

(a)_n = a(a-1)(a-2)...(a-n+1)

1385

(a)_n = a(a-1)(a-2)...(a-n+1)

\end{equation}

1386

\end{equation}

1387

The generalization of Pochhammer's symbol is given by equation \ref{poch}

1388

The generalization of Pochhammer's symbol is given by equation \ref{poch}

\begin{equation}

1389

\begin{equation}

\label{poch}

1390

\label{poch}

poch(a,x)= (a)_x = {\Gamma(a+x) \over \Gamma(a)}

1391

poch(a,x)= (a)_x = {\Gamma(a+x) \over \Gamma(a)}

\end{equation}

1392

\end{equation}

1393

Specific interfaces : \verb'poch,dpoch'

1394

Specific interfaces : \verb'poch,dpoch'

1395

1396

\subsubsection{poch1}

1397

\subsubsection{poch1}

\begin{verbatim}

1398

\begin{verbatim}

poch1(a,x)

1399

poch1(a,x)

\end{verbatim}

1400

\end{verbatim}

\begin{itemize}

1401

\begin{itemize}

\item x is a real

1402

\item x is a real

\item a is a real

1403

\item a is a real

\end{itemize}

1404

\end{itemize}

This function is defined by equation \ref{poch1}. It is usefull for certains situations, especially when x is small.

1405

This function is defined by equation \ref{poch1}. It is usefull for certains situations, especially when x is small.

1406

\begin{equation}

1407

\begin{equation}

\label{poch1}

1408

\label{poch1}

poch1(a,x)={{(a)_x-1} \over x}

1409

poch1(a,x)={{(a)_x-1} \over x}

\end{equation}

1410

\end{equation}

1411

Specific interfaces : \verb'poch1,dpoch1'

1412

Specific interfaces : \verb'poch1,dpoch1'

1413

\subsubsection{beta}

1414

\subsubsection{beta}

\begin{verbatim}

1415

\begin{verbatim}

beta(a,b)

1416

beta(a,b)

\end{verbatim}

1417

\end{verbatim}

\begin{itemize}

1418

\begin{itemize}

\item a,b are real positive or complex

1419

\item a,b are real positive or complex

\end{itemize}

1420

\end{itemize}

This function evaluates $\beta$ function defined by equation \ref{beta}.

1421

This function evaluates $\beta$ function defined by equation \ref{beta}.

\begin{equation}

1422

\begin{equation}

\label{beta}

1423

\label{beta}

beta(a,b)=\beta(a,b)={ {\Gamma(a) \Gamma(b)} \over \Gamma(a+b) }

1424

beta(a,b)=\beta(a,b)={ {\Gamma(a) \Gamma(b)} \over \Gamma(a+b) }

\end{equation}

1425

\end{equation}

1426

Specific interfaces : \verb'beta,dbeta,cbeta,zbeta'

1427

Specific interfaces : \verb'beta,dbeta,cbeta,zbeta'

1428

1429

\subsubsection{albeta}

1430

\subsubsection{albeta}

\begin{verbatim}

1431

\begin{verbatim}

albeta(a,b)

1432

albeta(a,b)

\end{verbatim}

1433

\end{verbatim}

\begin{itemize}

1434

\begin{itemize}

\item a,b are real positive or complex

1435

\item a,b are real positive or complex

\end{itemize}

1436

\end{itemize}

This function evaluates the natural logarithm of beta function : $ln(\beta(a,b))$

1437

This function evaluates the natural logarithm of beta function : $ln(\beta(a,b))$

1438

Specific interfaces : \verb'albeta,dlbeta,clbeta,zlbeta'

1439

Specific interfaces : \verb'albeta,dlbeta,clbeta,zlbeta'

1440

\subsubsection{betai}

1441

\subsubsection{betai}

\begin{verbatim}

1442

\begin{verbatim}

betai(x,pin,qin)

1443

betai(x,pin,qin)

\end{verbatim}

1444

\end{verbatim}

\begin{itemize}

1445

\begin{itemize}

\item x is a real in [0,1]

1446

\item x is a real in [0,1]

\item pin and qin are strictly positive real

1447

\item pin and qin are strictly positive real

\end{itemize}

1448

\end{itemize}

This function evaluates the incomplete beta function ratio, that is the probability that a random variable from a beta distribution having parameters pin and qin will be less than or equal to x. It is defined by equation \ref{betai}.

1449

This function evaluates the incomplete beta function ratio, that is the probability that a random variable from a beta distribution having parameters pin and qin will be less than or equal to x. It is defined by equation \ref{betai}.

1450

\begin{equation}

1451

\begin{equation}

\label{betai}

1452

\label{betai}

betai(x,pin,qin)=I_x(pin,qin)={1 \over \beta(pin,qin)} \int _0 ^x t^{pin-1}(1-t)^{qin-1}dt

1453

betai(x,pin,qin)=I_x(pin,qin)={1 \over \beta(pin,qin)} \int _0 ^x t^{pin-1}(1-t)^{qin-1}dt

\end{equation}

1454

\end{equation}

1455

Specific interfaces : \verb'betai,dbetai'

1456

Specific interfaces : \verb'betai,dbetai'

1457

\subsection{Error function and related}

1458

\subsection{Error function and related}

\subsubsection{erf}

1459

\subsubsection{erf}

\begin{verbatim}

1460

\begin{verbatim}

erf(x)

1461

erf(x)

\end{verbatim}

1462

\end{verbatim}

\begin{itemize}

1463

\begin{itemize}

\item x is a real

1464

\item x is a real

\end{itemize}

1465

\end{itemize}

This function evaluates the error function defined by equation \ref{erf}.

1466

This function evaluates the error function defined by equation \ref{erf}.

\begin{equation}

1467

\begin{equation}

\label{erf}

1468

\label{erf}

erf(x)={2\over \sqrt{ \pi}} \int _0 ^x e^{-t^2}dt

1469

erf(x)={2\over \sqrt{ \pi}} \int _0 ^x e^{-t^2}dt

\end{equation}

1470

\end{equation}

1471

Specific interfaces : \verb'erf,derf'

1472

Specific interfaces : \verb'erf,derf'

1473

\subsubsection{erfc}

1474

\subsubsection{erfc}

\begin{verbatim}

1475

\begin{verbatim}

erfc(x)

1476

erfc(x)

\end{verbatim}

1477

\end{verbatim}

\begin{itemize}

1478

\begin{itemize}

\item x is a real

1479

\item x is a real

\end{itemize}

1480

\end{itemize}

This function evaluates the complimentary error function defined by equation \ref{erfc}.

1481

This function evaluates the complimentary error function defined by equation \ref{erfc}.

\begin{equation}

1482

\begin{equation}

\label{erfc}

1483

\label{erfc}

erfc(x)={2\over \sqrt{ \pi}} \int _x ^\infty e^{-t^2}dt

1484

erfc(x)={2\over \sqrt{ \pi}} \int _x ^\infty e^{-t^2}dt

\end{equation}

1485

\end{equation}

1486

Specific interfaces : \verb'erfc,derfc'

1487

Specific interfaces : \verb'erfc,derfc'

1488

1489

\subsubsection{daws}

1490

\subsubsection{daws}

\begin{verbatim}

1491

\begin{verbatim}

daws(x)

1492

daws(x)

\end{verbatim}

1493

\end{verbatim}

\begin{itemize}

1494

\begin{itemize}

\item x is a real

1495

\item x is a real

\end{itemize}

1496

\end{itemize}

This function evaluates Dawson's function defined by equation \ref{daws}.

1497

This function evaluates Dawson's function defined by equation \ref{daws}.

\begin{equation}

1498

\begin{equation}

\label{daws}

1499

\label{daws}

daws(x)=e^{-x^2} \int _0 ^x e^{t^2}dt

1500

daws(x)=e^{-x^2} \int _0 ^x e^{t^2}dt

\end{equation}

1501

\end{equation}

1502

Specific interfaces : \verb'daws,ddaws'

1503

Specific interfaces : \verb'daws,ddaws'

1504

\subsection{Bessel functions and related}

1505

\subsection{Bessel functions and related}

\subsubsection{bsj0}

1506

\subsubsection{bsj0}

\begin{verbatim}

1507

\begin{verbatim}

bsj0(x)

1508

bsj0(x)

\end{verbatim}

1509

\end{verbatim}

\begin{itemize}

1510

\begin{itemize}

\item x is a real

1511

\item x is a real

\end{itemize}

1512

\end{itemize}

This function evaluates Bessel function of the first kind of order 0 defined by equation \ref{bsj0}.

1513

This function evaluates Bessel function of the first kind of order 0 defined by equation \ref{bsj0}.

\begin{equation}

1514

\begin{equation}

\label{bsj0}

1515

\label{bsj0}

bsj0(x)=J_0(x)= {1 \over \pi} \int _0 ^\pi cos(x sin(\theta)) d\theta

1516

bsj0(x)=J_0(x)= {1 \over \pi} \int _0 ^\pi cos(x sin(\theta)) d\theta

\end{equation}

1517

\end{equation}

1518

Specific interfaces : \verb'besj0,dbesj0'

1519

Specific interfaces : \verb'besj0,dbesj0'

1520

GITLAB

git-svn-id: https://lxsd.femto-st.fr/svn/fvn@35 b657c933-2333-4658-acf2-d3c7c2708721