Commit 59ae88e0693ec0d9934298af451f5e17c5e713d0

Authored by daniau
1 parent 02a8cc89cd

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

Showing 8 changed files with 65 additions and 24 deletions Inline Diff

1 1
include $(BTREE)/Make.inc 2 2 include $(BTREE)/Make.inc
3 3
4 4
objects = fvnlib.o 5 5 objects = fvnlib.o
library = libfvn.a 6 6 library = libfvn$(libext)
7 7
all: umfpack fnlib $(library) 8 8 all: umfpack fnlib $(library)
9 9
clean: 10 10 clean:
rm -f {*.o,*.oo,*.ipo,*.a,*.mod} 11 11 rm -f {*.o,*.oo,*.ipo,*.a,*.mod}
( cd fvn_sparse ; rm -f {*.o,*.oo,*.ipo,*.a,*.mod} ) 12 12 ( cd fvn_sparse ; rm -f {*.o,*.oo,*.ipo,*.a,*.mod} )
( cd fvn_sparse/AMD ; make clean ) 13 13 ( cd fvn_sparse/AMD ; make clean )
( cd fvn_sparse/UMFPACK ; make clean ) 14 14 ( cd fvn_sparse/UMFPACK ; make clean )
( cd fvn_fnlib ; make clean ) 15 15 ( cd fvn_fnlib ; make clean )
( rm -f fvn_sparse/AMD/Lib/libamd.a ) 16 16 ( rm -f fvn_sparse/AMD/Lib/libamd.a )
( rm -f fvn_sparse/UMFPACK/Lib/libumfpack.a ) 17 17 ( rm -f fvn_sparse/UMFPACK/Lib/libumfpack.a )
18 18
install: 19 19 install:
cp fvn.mod $(BTREE)/modules 20 20 cp fvn.mod $(BTREE)/modules
cp libfvn.a $(BTREE)/lib 21 21 cp $(library) $(BTREE)/lib
( cd fvn_fnlib ; make install ) 22 22 ( cd fvn_fnlib ; make install )
( cp fvn_sparse/UMFPACK/Lib/libumfpack.a $(BTREE)/lib ) 23 23 ( cp fvn_sparse/UMFPACK/Lib/libumfpack.a $(BTREE)/lib/libumfpack.lib )
( cp fvn_sparse/AMD/Lib/libamd.a $(BTREE)/lib ) 24 24 ( cp fvn_sparse/AMD/Lib/libamd.a $(BTREE)/lib/libamd.lib )
25 25
%.o: %.f90 26 26 %.o: %.f90
$(F95) $(F95FLAGS) -c $< 27 27 $(F95) $(F95FLAGS) -c $< -o $@
28 28
$(objects):fvnlib.f90 fvn_quadpack/dqk15_2d_inner.f fvn_quadpack/dqk31_2d_outer.f \ 29 29 $(objects):fvnlib.f90 fvn_quadpack/dqk15_2d_inner.f fvn_quadpack/dqk31_2d_outer.f \
fvn_quadpack/d1mach.f fvn_quadpack/dqk31_2d_inner.f fvn_quadpack/dqage.f \ 30 30 fvn_quadpack/d1mach.f fvn_quadpack/dqk31_2d_inner.f fvn_quadpack/dqage.f \
fvn_quadpack/dqk15.f fvn_quadpack/dqk21.f fvn_quadpack/dqk31.f fvn_quadpack/dqk41.f \ 31 31 fvn_quadpack/dqk15.f fvn_quadpack/dqk21.f fvn_quadpack/dqk31.f fvn_quadpack/dqk41.f \
fvn_quadpack/dqk51.f fvn_quadpack/dqk61.f fvn_quadpack/dqk41_2d_outer.f \ 32 32 fvn_quadpack/dqk51.f fvn_quadpack/dqk61.f fvn_quadpack/dqk41_2d_outer.f \
fvn_quadpack/dqk41_2d_inner.f fvn_quadpack/dqag_2d_outer.f fvn_quadpack/dqag_2d_inner.f \ 33 33 fvn_quadpack/dqk41_2d_inner.f fvn_quadpack/dqag_2d_outer.f fvn_quadpack/dqag_2d_inner.f \
fvn_quadpack/dqpsrt.f fvn_quadpack/dqag.f fvn_quadpack/dqage_2d_outer.f \ 34 34 fvn_quadpack/dqpsrt.f fvn_quadpack/dqag.f fvn_quadpack/dqage_2d_outer.f \
fvn_quadpack/dqage_2d_inner.f fvn_quadpack/dqk51_2d_outer.f fvn_quadpack/dqk51_2d_inner.f \ 35 35 fvn_quadpack/dqage_2d_inner.f fvn_quadpack/dqk51_2d_outer.f fvn_quadpack/dqk51_2d_inner.f \
fvn_quadpack/dqk61_2d_outer.f fvn_quadpack/dqk21_2d_outer.f fvn_quadpack/dqk61_2d_inner.f \ 36 36 fvn_quadpack/dqk61_2d_outer.f fvn_quadpack/dqk21_2d_outer.f fvn_quadpack/dqk61_2d_inner.f \
fvn_quadpack/dqk21_2d_inner.f fvn_quadpack/dqk15_2d_outer.f 37 37 fvn_quadpack/dqk21_2d_inner.f fvn_quadpack/dqk15_2d_outer.f
38 38
$(library): $(objects) 39 39 $(library): $(objects)
$(AR) rcu $@ $(objects) 40 40 $(AR) rcu $@ $(objects)
$(AR) ru $@ fvn_sparse/umfpack_wrapper.o 41 41 $(AR) ru $@ fvn_sparse/umfpack_wrapper.o
$(AR) s $@ 42 42 $(RANLIB) $@
43 43
umfpack: 44 44 umfpack:
( cd fvn_sparse/UMFPACK ; make ) 45 45 ( cd fvn_sparse/UMFPACK ; make )
1 1
include $(BTREE)/Make.inc 2 2 include $(BTREE)/Make.inc
3 3
library = libfvn_fnlib.a 4 4 library = libfvn_fnlib$(libext)
5 5
objects = acosh.o aide.o aid.o aie.o \ 6 6 objects = acosh.o aide.o aid.o aie.o \
ai.o albeta.o algams.o ali.o \ 7 7 ai.o albeta.o algams.o ali.o \
alngam.o alnrel.o asinh.o atanh.o \ 8 8 alngam.o alnrel.o asinh.o atanh.o \
besi0e.o besi0.o besi1e.o besi1.o \ 9 9 besi0e.o besi0.o besi1e.o besi1.o \
besj0.o besj1.o besk0e.o besk0.o \ 10 10 besj0.o besj1.o besk0e.o besk0.o \
besk1e.o besk1.o beskes.o besks.o \ 11 11 besk1e.o besk1.o beskes.o besks.o \
besy0.o besy1.o betai.o beta.o \ 12 12 besy0.o besy1.o betai.o beta.o \
bide.o bid.o bie.o binom.o \ 13 13 bide.o bid.o bie.o binom.o \
bi.o c0lgmc.o c8lgmc.o c9lgmc.o \ 14 14 bi.o c0lgmc.o c8lgmc.o c9lgmc.o \
c9ln2r.o cacosh.o cacos.o carg.o \ 15 15 c9ln2r.o cacosh.o cacos.o carg.o \
casinh.o casin.o catan2.o catanh.o \ 16 16 casinh.o casin.o catan2.o catanh.o \
catan.o cbeta.o cbrt.o ccbrt.o \ 17 17 catan.o cbeta.o cbrt.o ccbrt.o \
ccosh.o ccot.o cexprl.o cgamma.o \ 18 18 ccosh.o ccot.o cexprl.o cgamma.o \
cgamr.o chi.o chu.o cinh.o \ 19 19 cgamr.o chi.o chu.o cinh.o \
cin.o ci.o clbeta.o clngam.o \ 20 20 cin.o ci.o clbeta.o clngam.o \
clnrel.o clog10.o comp1.o comp2.o \ 21 21 clnrel.o clog10.o comp1.o comp2.o \
comp3.o cosdg.o cot.o cpsi.o \ 22 22 comp3.o cosdg.o cot.o cpsi.o \
csevl.o csinh.o ctanh.o ctan.o \ 23 23 csevl.o csinh.o ctanh.o ctan.o \
d1mach.o d9admp.o d9aimp.o d9atn1.o \ 24 24 d1mach.o d9admp.o d9aimp.o d9atn1.o \
d9b0mp.o d9b1mp.o d9chm.o d9chu.o \ 25 25 d9b0mp.o d9b1mp.o d9chm.o d9chu.o \
d9gaml.o d9gmic.o d9gmit.o d9knus.o \ 26 26 d9gaml.o d9gmic.o d9gmit.o d9knus.o \
d9lgic.o d9lgit.o d9lgmc.o d9ln2r.o \ 27 27 d9lgic.o d9lgit.o d9lgmc.o d9ln2r.o \
d9pak.o d9sifg.o d9upak.o dacosh.o \ 28 28 d9pak.o d9sifg.o d9upak.o dacosh.o \
daide.o daid.o daie.o dai.o \ 29 29 daide.o daid.o daie.o dai.o \
dasinh.o datanh.o daws.o dbesi0.o \ 30 30 dasinh.o datanh.o daws.o dbesi0.o \
dbesi1.o dbesj0.o dbesj1.o dbesk0.o \ 31 31 dbesi1.o dbesj0.o dbesj1.o dbesk0.o \
dbesk1.o dbesks.o dbesy0.o dbesy1.o \ 32 32 dbesk1.o dbesks.o dbesy0.o dbesy1.o \
dbetai.o dbeta.o dbide.o dbid.o \ 33 33 dbetai.o dbeta.o dbide.o dbid.o \
dbie.o dbinom.o dbi.o dbsi0e.o \ 34 34 dbie.o dbinom.o dbi.o dbsi0e.o \
dbsi1e.o dbsk0e.o dbsk1e.o dbskes.o \ 35 35 dbsi1e.o dbsk0e.o dbsk1e.o dbskes.o \
dcbrt.o dchi.o dchu.o dcinh.o \ 36 36 dcbrt.o dchi.o dchu.o dcinh.o \
dcin.o dci.o dcosdg.o dcot.o \ 37 37 dcin.o dci.o dcosdg.o dcot.o \
dcsevl.o ddaws.o de1.o dei.o \ 38 38 dcsevl.o ddaws.o de1.o dei.o \
derfc.o derf.o dexprl.o dfac.o \ 39 39 derfc.o derf.o dexprl.o dfac.o \
dgamic.o dgami.o dgamit.o dgamma.o \ 40 40 dgamic.o dgami.o dgamit.o dgamma.o \
dgamr.o dlbeta.o dlgams.o dli.o \ 41 41 dgamr.o dlbeta.o dlgams.o dli.o \
dlngam.o dlnrel.o dpoch1.o dpoch.o \ 42 42 dlngam.o dlnrel.o dpoch1.o dpoch.o \
dpsi.o dshi.o dsindg.o dsi.o \ 43 43 dpsi.o dshi.o dsindg.o dsi.o \
dspenc.o e1.o e9rint.o ei.o \ 44 44 dspenc.o e1.o e9rint.o ei.o \
entsrc.o eprint.o erfc.o erf.o \ 45 45 entsrc.o eprint.o erfc.o erf.o \
erroff.o exprel.o fac.o fdump.o \ 46 46 erroff.o exprel.o fac.o fdump.o \
fvn_fnlib.o gamic.o gami.o gamit.o \ 47 47 fvn_fnlib.o gamic.o gami.o gamit.o \
gamma.o gamr.o i1mach.o i8save.o \ 48 48 gamma.o gamr.o i1mach.o i8save.o \
initds.o inits.o nerror.o poch1.o \ 49 49 initds.o inits.o nerror.o poch1.o \
poch.o psi.o r1mach.o r9admp.o \ 50 50 poch.o psi.o r1mach.o r9admp.o \
r9aimp.o r9atn1.o r9chm.o r9chu.o \ 51 51 r9aimp.o r9atn1.o r9chm.o r9chu.o \
r9gaml.o r9gmic.o r9gmit.o r9knus.o \ 52 52 r9gaml.o r9gmic.o r9gmit.o r9knus.o \
r9lgic.o r9lgit.o r9lgmc.o r9ln2r.o \ 53 53 r9lgic.o r9lgit.o r9lgmc.o r9ln2r.o \
r9pak.o r9sifg.o r9upak.o randgs.o \ 54 54 r9pak.o r9sifg.o r9upak.o randgs.o \
rand.o random.o ranf.o retsrc.o \ 55 55 rand.o random.o ranf.o retsrc.o \
s88fmt.o s9comp.o seterr.o seteru.o \ 56 56 s88fmt.o s9comp.o seterr.o seteru.o \
shi.o sindg.o si.o spenc.o \ 57 57 shi.o sindg.o si.o spenc.o \
z0lgmc.o z8lgmc.o z9lgmc.o z9ln2r.o \ 58 58 z0lgmc.o z8lgmc.o z9lgmc.o z9ln2r.o \
zacosh.o zacos.o zarg.o zasinh.o \ 59 59 zacosh.o zacos.o zarg.o zasinh.o \
zasin.o zatan2.o zatanh.o zatan.o \ 60 60 zasin.o zatan2.o zatanh.o zatan.o \
zbeta.o zcbrt.o zcosh.o zcot.o \ 61 61 zbeta.o zcbrt.o zcosh.o zcot.o \
zexprl.o zgamma.o zgamr.o zlbeta.o \ 62 62 zexprl.o zgamma.o zgamr.o zlbeta.o \
zlngam.o zlnrel.o zlog10.o zpsi.o \ 63 63 zlngam.o zlnrel.o zlog10.o zpsi.o \
zsinh.o ztanh.o ztan.o besyn.o \ 64 64 zsinh.o ztanh.o ztan.o besyn.o \
besjn.o dbesyn.o dbesjn.o beskn.o \ 65 65 besjn.o dbesyn.o dbesjn.o beskn.o \
besin.o dbeskn.o dbesin.o 66 66 besin.o dbeskn.o dbesin.o
67 67
lib:$(library) 68 68 lib:$(library)
69 69
$(library): $(objects) 70 70 $(library): $(objects)
$(AR) rcu $@ $(objects) 71 71 $(AR) rcu $@ $(objects)
$(AR) s $@ 72 72 $(RANLIB) $@
73 73
install: 74 74 install:
cp fvn_fnlib.mod $(BTREE)/modules 75 75 cp fvn_fnlib.mod $(BTREE)/modules
cp libfvn_fnlib.a $(BTREE)/lib 76 76 cp $(library) $(BTREE)/lib
77 77
clean: 78 78 clean:
rm -f {*.o,*.oo,*.ipo,*.a,*.mod} 79 79 rm -f {*.o,*.oo,*.ipo,*.a,*.mod}
80 80
%.o: %.f90 81 81 %.o: %.f90
function alngam (x) 1 1 function alngam (x)
c august 1980 edition. w. fullerton, c3, los alamos scientific lab. 2 2 c august 1980 edition. w. fullerton, c3, los alamos scientific lab.
c 3 3 c
external gamma, r1mach, r9lgmc 4 4 external gamma, r1mach, r9lgmc
data sq2pil / 0.9189385332 0467274e0/ 5 5 data sq2pil / 0.9189385332 0467274e0/
c sq2pil = alog(sqrt(2.*pi)), sqpi2l = alog (sqrt(pi/2.)) 6 6 c sq2pil = alog(sqrt(2.*pi)), sqpi2l = alog (sqrt(pi/2.))
data sqpi2l / 0.2257913526 4472743e0/ 7 7 data sqpi2l / 0.2257913526 4472743e0/
data pi / 3.1415926535 8979324e0/ 8 8 data pi / 3.1415926535 8979324e0/
c 9 9 c
data xmax, dxrel / 0., 0. / 10 10 data xmax, dxrel / 0., 0. /
c 11 11 c
if (xmax.ne.0.) go to 10 12 12 if (xmax.ne.0.) go to 10
xmax = r1mach(2)/alog(r1mach(2)) 13 13 xmax = r1mach(2)/alog(r1mach(2))
dxrel = sqrt (r1mach(4)) 14 14 dxrel = sqrt (r1mach(4))
c 15 15 c
10 y = abs(x) 16 16 10 y = abs(x)
if (y.gt.10.0) go to 20 17 17 if (y.gt.10.0) go to 20
c 18 18 c
c alog (abs (gamma(x))) for abs(x) .le. 10.0 19 19 c alog (abs (gamma(x))) for abs(x) .le. 10.0
c 20 20 c
alngam = alog (abs (gamma(x))) 21 21 alngam = alog (abs (gamma(x)))
return 22 22 return
c 23 23 c
c alog (abs (gamma(x))) for abs(x) .gt. 10.0 24 24 c alog (abs (gamma(x))) for abs(x) .gt. 10.0
c 25 25 c
20 if (y.gt.xmax) call seteru ( 26 26 20 if (y.gt.xmax) call seteru (
1 38halngam abs(x) so big alngam overflows, 38, 2, 2) 27 27 1 38halngam abs(x) so big alngam overflows, 38, 2, 2)
c 28 28 c
if (x.gt.0.) alngam = sq2pil + (x-0.5)*alog(x) - x + r9lgmc(y) 29 29 if (x.gt.0.) alngam = sq2pil + (x-0.5)*alog(x) - x + r9lgmc(y)
if (x.gt.0.) return 30 30 if (x.gt.0.) return
c 31 31 c
sinpiy = abs (sin(pi*y)) 32 32 sinpiy = abs (sin(pi*y))
if (sinpiy.eq.0.) call seteru (31halngam x is a negative integer, 33 33 if (sinpiy.eq.0.) call seteru (31halngam x is a negative integer,
1 31, 3, 2) 34 34 1 31, 3, 2)
c 35 35 c
alngam = sqpi2l + (x-0.5)*alog(y) - x - alog(sinpiy) - r9lgmc(y) 36 36 alngam = sqpi2l + (x-0.5)*alog(y) - x - alog(sinpiy) - r9lgmc(y)
c 37 37 c
if (abs((x-aint(x-0.5))*alngam/x).lt.dxrel) call seteru ( 38 38 if (abs((x-aint(x-0.5))*alngam/x).lt.dxrel) call seteru (
1'alngam answer lt half precision because x too near negative ', 39 39 160halngam answer lt half precision because x too near negative,
2 68, 1, 1) 40 40 2 68, 1, 1)
return 41 41 return
c 42 42 c
end 43 43 end
44 44
fvn_sparse/UFconfig/UFconfig.mk
# This is a fvn specific UFconfig.mk file 1 1 # This is a fvn specific UFconfig.mk file
2 2
3 3
include $(BTREE)/Make.inc 4 4 include $(BTREE)/Make.inc
5 5
# These are defined in Make.inc but could be overriden here 6 6 # These are defined in Make.inc but could be overriden here
#CC = gcc 7 7 #CC = gcc
#CFLAGS = -O2 -m64 -fpic 8 8 #CFLAGS = -O2 -m64 -fpic
9 9
# ranlib, and ar, for generating libraries 10 10 # ranlib, and ar, for generating libraries
RANLIB = ranlib 11 11 RANLIB = ranlib
AR = ar cr 12 12 AR = ar cr
13 13
# delete and rename a file 14 14 # delete and rename a file
RM = rm -f 15 15 RM = rm -f
MV = mv -f 16 16 MV = mv -f
17 17
# Fortran compiler (not normally required) 18 18 # Fortran compiler (not normally required)
F77 = $(F95) 19 19 F77 = $(F95)
F77FLAGS = $(F95FLAGS) 20 20 F77FLAGS = $(F95FLAGS)
F77LIB = 21 21 F77LIB =
22 22
# C and Fortran libraries 23 23 # C and Fortran libraries
# No use to repeat BLAS here as it should be defined in Make.inc 24 24 # No use to repeat BLAS here as it should be defined in Make.inc
# but in case we are using acml we include it here 25 25 # but in case we are using acml we include it here
LIB = -lm $(ACML) $(EXTRALIBS) 26 26 LIB = $(ACML) $(EXTRALIBS)
27 27
XERBLA = 28 28 XERBLA =
29 29
UMFPACK_CONFIG = 30 30 UMFPACK_CONFIG =
31 31
#------------------------------------------------------------------------------ 32 32 #------------------------------------------------------------------------------
# Linux 33 33 # Linux
#------------------------------------------------------------------------------ 34 34 #------------------------------------------------------------------------------
35 35
# Using default compilers: 36 36 # Using default compilers:
# CC = gcc 37 37 # CC = gcc
# CFLAGS = -O3 38 38 # CFLAGS = -O3
39 39
# alternatives: 40 40 # alternatives:
# CFLAGS = -g -fexceptions \ 41 41 # CFLAGS = -g -fexceptions \
-Wall -W -Wshadow -Wmissing-prototypes -Wstrict-prototypes \ 42 42 -Wall -W -Wshadow -Wmissing-prototypes -Wstrict-prototypes \
-Wredundant-decls -Wnested-externs -Wdisabled-optimization -ansi 43 43 -Wredundant-decls -Wnested-externs -Wdisabled-optimization -ansi
# CFLAGS = -O3 -fexceptions \ 44 44 # CFLAGS = -O3 -fexceptions \
-Wall -W -Werror -Wshadow -Wmissing-prototypes -Wstrict-prototypes \ 45 45 -Wall -W -Werror -Wshadow -Wmissing-prototypes -Wstrict-prototypes \
-Wredundant-decls -Wnested-externs -Wdisabled-optimization -ansi 46 46 -Wredundant-decls -Wnested-externs -Wdisabled-optimization -ansi
#CFLAGS = -O3 -fexceptions -D_FILE_OFFSET_BITS=64 -D_LARGEFILE64_SOURCE 47 47 #CFLAGS = -O3 -fexceptions -D_FILE_OFFSET_BITS=64 -D_LARGEFILE64_SOURCE
# CFLAGS = -O3 48 48 # CFLAGS = -O3
49 49
# consider: 50 50 # consider:
# -fforce-addr -fmove-all-movables -freduce-all-givs -ftsp-ordering 51 51 # -fforce-addr -fmove-all-movables -freduce-all-givs -ftsp-ordering
# -frename-registers -ffast-math -funroll-loops 52 52 # -frename-registers -ffast-math -funroll-loops
53 53
# Using the Goto BLAS: 54 54 # Using the Goto BLAS:
# BLAS = -lgoto -lfrtbegin -lg2c $(XERBLA) -lpthread 55 55 # BLAS = -lgoto -lfrtbegin -lg2c $(XERBLA) -lpthread
56 56
# Using Intel's icc and ifort compilers: 57 57 # Using Intel's icc and ifort compilers:
# (does not work for mexFunctions unless you add a mexopts.sh file) 58 58 # (does not work for mexFunctions unless you add a mexopts.sh file)
# F77 = ifort 59 59 # F77 = ifort
# CC = icc 60 60 # CC = icc
# CFLAGS = -O3 -xN -vec_report=0 61 61 # CFLAGS = -O3 -xN -vec_report=0
# CFLAGS = -g 62 62 # CFLAGS = -g
fvn_sparse/UMFPACK/Source/cholmod_blas.h
/* ========================================================================== */ 1 1 /* ========================================================================== */
/* === Include/cholmod_blas.h =============================================== */ 2 2 /* === Include/cholmod_blas.h =============================================== */
/* ========================================================================== */ 3 3 /* ========================================================================== */
4 4
/* ----------------------------------------------------------------------------- 5 5 /* -----------------------------------------------------------------------------
* CHOLMOD/Include/cholmod_blas.h. 6 6 * CHOLMOD/Include/cholmod_blas.h.
* Copyright (C) Univ. of Florida. Author: Timothy A. Davis 7 7 * Copyright (C) Univ. of Florida. Author: Timothy A. Davis
* CHOLMOD/Include/cholmod_blas.h is licensed under Version 2.1 of the GNU 8 8 * CHOLMOD/Include/cholmod_blas.h is licensed under Version 2.1 of the GNU
* Lesser General Public License. See lesser.txt for a text of the license. 9 9 * Lesser General Public License. See lesser.txt for a text of the license.
* CHOLMOD is also available under other licenses; contact authors for details. 10 10 * CHOLMOD is also available under other licenses; contact authors for details.
* http://www.cise.ufl.edu/research/sparse 11 11 * http://www.cise.ufl.edu/research/sparse
* -------------------------------------------------------------------------- */ 12 12 * -------------------------------------------------------------------------- */
13 13
/* This does not need to be included in the user's program. */ 14 14 /* This does not need to be included in the user's program. */
15 15
#ifndef CHOLMOD_BLAS_H 16 16 #ifndef CHOLMOD_BLAS_H
#define CHOLMOD_BLAS_H 17 17 #define CHOLMOD_BLAS_H
18 18
/* ========================================================================== */ 19 19 /* ========================================================================== */
/* === Architecture ========================================================= */ 20 20 /* === Architecture ========================================================= */
/* ========================================================================== */ 21 21 /* ========================================================================== */
22 22
#if defined (__sun) || defined (MSOL2) || defined (ARCH_SOL2) 23 23 #if defined (__sun) || defined (MSOL2) || defined (ARCH_SOL2)
#define CHOLMOD_SOL2 24 24 #define CHOLMOD_SOL2
#define CHOLMOD_ARCHITECTURE "Sun Solaris" 25 25 #define CHOLMOD_ARCHITECTURE "Sun Solaris"
26 26
#elif defined (__sgi) || defined (MSGI) || defined (ARCH_SGI) 27 27 #elif defined (__sgi) || defined (MSGI) || defined (ARCH_SGI)
#define CHOLMOD_SGI 28 28 #define CHOLMOD_SGI
#define CHOLMOD_ARCHITECTURE "SGI Irix" 29 29 #define CHOLMOD_ARCHITECTURE "SGI Irix"
30 30
#elif defined (__linux) || defined (MGLNX86) || defined (ARCH_GLNX86) 31 31 #elif defined (__linux) || defined (MGLNX86) || defined (ARCH_GLNX86)
#define CHOLMOD_LINUX 32 32 #define CHOLMOD_LINUX
#define CHOLMOD_ARCHITECTURE "Linux" 33 33 #define CHOLMOD_ARCHITECTURE "Linux"
34 34
#elif defined (_AIX) || defined (MIBM_RS) || defined (ARCH_IBM_RS) 35 35 #elif defined (_AIX) || defined (MIBM_RS) || defined (ARCH_IBM_RS)
#define CHOLMOD_AIX 36 36 #define CHOLMOD_AIX
#define CHOLMOD_ARCHITECTURE "IBM AIX" 37 37 #define CHOLMOD_ARCHITECTURE "IBM AIX"
#define BLAS_NO_UNDERSCORE 38 38 #define BLAS_NO_UNDERSCORE
39 39
#elif defined (__alpha) || defined (MALPHA) || defined (ARCH_ALPHA) 40 40 #elif defined (__alpha) || defined (MALPHA) || defined (ARCH_ALPHA)
#define CHOLMOD_ALPHA 41 41 #define CHOLMOD_ALPHA
#define CHOLMOD_ARCHITECTURE "Compaq Alpha" 42 42 #define CHOLMOD_ARCHITECTURE "Compaq Alpha"
43 43
#elif defined (_WIN32) || defined (WIN32) || defined (_WIN64) || defined (WIN64) 44 44 #elif defined (_WIN32) || defined (WIN32) || defined (_WIN64) || defined (WIN64)
#if defined (__MINGW32__) || defined (__MINGW32__) 45 45 #if defined (__MINGW32__) || defined (__MINGW32__)
#define CHOLMOD_MINGW 46 46 #define CHOLMOD_MINGW
#elif defined (__CYGWIN32__) || defined (__CYGWIN32__) 47 47 #elif defined (__CYGWIN32__) || defined (__CYGWIN32__)
#define CHOLMOD_CYGWIN 48 48 #define CHOLMOD_CYGWIN
#else 49 49 #else
#define CHOLMOD_WINDOWS 50 50 #define CHOLMOD_WINDOWS
#define BLAS_NO_UNDERSCORE 51 51 #define BLAS_NO_UNDERSCORE
#endif 52 52 #endif
#define CHOLMOD_ARCHITECTURE "Microsoft Windows" 53 53 #define CHOLMOD_ARCHITECTURE "Microsoft Windows"
54 54
#elif defined (__hppa) || defined (__hpux) || defined (MHPUX) || defined (ARCH_HPUX) 55 55 #elif defined (__hppa) || defined (__hpux) || defined (MHPUX) || defined (ARCH_HPUX)
#define CHOLMOD_HP 56 56 #define CHOLMOD_HP
#define CHOLMOD_ARCHITECTURE "HP Unix" 57 57 #define CHOLMOD_ARCHITECTURE "HP Unix"
#define BLAS_NO_UNDERSCORE 58 58 #define BLAS_NO_UNDERSCORE
59 59
#elif defined (__hp700) || defined (MHP700) || defined (ARCH_HP700) 60 60 #elif defined (__hp700) || defined (MHP700) || defined (ARCH_HP700)
#define CHOLMOD_HP 61 61 #define CHOLMOD_HP
#define CHOLMOD_ARCHITECTURE "HP 700 Unix" 62 62 #define CHOLMOD_ARCHITECTURE "HP 700 Unix"
#define BLAS_NO_UNDERSCORE 63 63 #define BLAS_NO_UNDERSCORE
64 64
#else 65 65 #else
/* If the architecture is unknown, and you call the BLAS, you may need to */ 66 66 /* If the architecture is unknown, and you call the BLAS, you may need to */
/* define BLAS_BY_VALUE, BLAS_NO_UNDERSCORE, and/or BLAS_CHAR_ARG yourself. */ 67 67 /* define BLAS_BY_VALUE, BLAS_NO_UNDERSCORE, and/or BLAS_CHAR_ARG yourself. */
#define CHOLMOD_ARCHITECTURE "unknown" 68 68 #define CHOLMOD_ARCHITECTURE "unknown"
#endif 69 69 #endif
70 70
71 71
/* ========================================================================== */ 72 72 /* ========================================================================== */
/* === BLAS and LAPACK names ================================================ */ 73 73 /* === BLAS and LAPACK names ================================================ */
/* ========================================================================== */ 74 74 /* ========================================================================== */
75 75
/* Prototypes for the various versions of the BLAS. */ 76 76 /* Prototypes for the various versions of the BLAS. */
77 77
/* Determine if the 64-bit Sun Performance BLAS is to be used */ 78 78 /* Determine if the 64-bit Sun Performance BLAS is to be used */
#if defined(CHOLMOD_SOL2) && !defined(NSUNPERF) && defined(LONG) && defined(LONGBLAS) 79 79 #if defined(CHOLMOD_SOL2) && !defined(NSUNPERF) && defined(LONG) && defined(LONGBLAS)
#define SUN64 80 80 #define SUN64
#endif 81 81 #endif
82 82
#ifdef SUN64 83 83 #ifdef SUN64
84 84
#define BLAS_DTRSV dtrsv_64_ 85 85 #define BLAS_DTRSV dtrsv_64_
#define BLAS_DGEMV dgemv_64_ 86 86 #define BLAS_DGEMV dgemv_64_
#define BLAS_DTRSM dtrsm_64_ 87 87 #define BLAS_DTRSM dtrsm_64_
#define BLAS_DGEMM dgemm_64_ 88 88 #define BLAS_DGEMM dgemm_64_
#define BLAS_DSYRK dsyrk_64_ 89 89 #define BLAS_DSYRK dsyrk_64_
#define BLAS_DGER dger_64_ 90 90 #define BLAS_DGER dger_64_
#define BLAS_DSCAL dscal_64_ 91 91 #define BLAS_DSCAL dscal_64_
#define LAPACK_DPOTRF dpotrf_64_ 92 92 #define LAPACK_DPOTRF dpotrf_64_
93 93
#define BLAS_ZTRSV ztrsv_64_ 94 94 #define BLAS_ZTRSV ztrsv_64_
#define BLAS_ZGEMV zgemv_64_ 95 95 #define BLAS_ZGEMV zgemv_64_
#define BLAS_ZTRSM ztrsm_64_ 96 96 #define BLAS_ZTRSM ztrsm_64_
#define BLAS_ZGEMM zgemm_64_ 97 97 #define BLAS_ZGEMM zgemm_64_
#define BLAS_ZHERK zherk_64_ 98 98 #define BLAS_ZHERK zherk_64_
#define BLAS_ZGER zgeru_64_ 99 99 #define BLAS_ZGER zgeru_64_
#define BLAS_ZSCAL zscal_64_ 100 100 #define BLAS_ZSCAL zscal_64_
#define LAPACK_ZPOTRF zpotrf_64_ 101 101 #define LAPACK_ZPOTRF zpotrf_64_
102 102
#elif defined (BLAS_NO_UNDERSCORE) 103 103 #elif defined (BLAS_NO_UNDERSCORE)
104 104
#define BLAS_DTRSV dtrsv 105 105 #define BLAS_DTRSV dtrsv
#define BLAS_DGEMV dgemv 106 106 #define BLAS_DGEMV dgemv
#define BLAS_DTRSM dtrsm 107 107 #define BLAS_DTRSM dtrsm
#define BLAS_DGEMM dgemm 108 108 #define BLAS_DGEMM dgemm
#define BLAS_DSYRK dsyrk 109 109 #define BLAS_DSYRK dsyrk
#define BLAS_DGER dger 110 110 #define BLAS_DGER dger
#define BLAS_DSCAL dscal 111 111 #define BLAS_DSCAL dscal
#define LAPACK_DPOTRF dpotrf 112 112 #define LAPACK_DPOTRF dpotrf
113 113
#define BLAS_ZTRSV ztrsv 114 114 #define BLAS_ZTRSV ztrsv
#define BLAS_ZGEMV zgemv 115 115 #define BLAS_ZGEMV zgemv
#define BLAS_ZTRSM ztrsm 116 116 #define BLAS_ZTRSM ztrsm
#define BLAS_ZGEMM zgemm 117 117 #define BLAS_ZGEMM zgemm
#define BLAS_ZHERK zherk 118 118 #define BLAS_ZHERK zherk
#define BLAS_ZGER zgeru 119 119 #define BLAS_ZGER zgeru
#define BLAS_ZSCAL zscal 120 120 #define BLAS_ZSCAL zscal
#define LAPACK_ZPOTRF zpotrf 121 121 #define LAPACK_ZPOTRF zpotrf
122 122
#else 123 123 #else
124 124
#define BLAS_DTRSV dtrsv_ 125 125 #define BLAS_DTRSV dtrsv_
#define BLAS_DGEMV dgemv_ 126 126 #define BLAS_DGEMV dgemv_
#define BLAS_DTRSM dtrsm_ 127 127 #define BLAS_DTRSM dtrsm_
#define BLAS_DGEMM dgemm_ 128 128 #define BLAS_DGEMM dgemm_
#define BLAS_DSYRK dsyrk_ 129 129 #define BLAS_DSYRK dsyrk_
#define BLAS_DGER dger_ 130 130 #define BLAS_DGER dger_
#define BLAS_DSCAL dscal_ 131 131 #define BLAS_DSCAL dscal_
#define LAPACK_DPOTRF dpotrf_ 132 132 #define LAPACK_DPOTRF dpotrf_
133 133
#define BLAS_ZTRSV ztrsv_ 134 134 #define BLAS_ZTRSV ztrsv_
#define BLAS_ZGEMV zgemv_ 135 135 #define BLAS_ZGEMV zgemv_
#define BLAS_ZTRSM ztrsm_ 136 136 #define BLAS_ZTRSM ztrsm_
#define BLAS_ZGEMM zgemm_ 137 137 #define BLAS_ZGEMM zgemm_
#define BLAS_ZHERK zherk_ 138 138 #define BLAS_ZHERK zherk_
#define BLAS_ZGER zgeru_ 139 139 #define BLAS_ZGER zgeru_
#define BLAS_ZSCAL zscal_ 140 140 #define BLAS_ZSCAL zscal_
#define LAPACK_ZPOTRF zpotrf_ 141 141 #define LAPACK_ZPOTRF zpotrf_
142 142
#endif 143 143 #endif
144 144
145 #ifdef PGIW32
146
147 #undef BLAS_DTRSV
148 #undef BLAS_DGEMV
149 #undef BLAS_DTRSM
150 #undef BLAS_DGEMM
151 #undef BLAS_DSYRK
152 #undef BLAS_DGER
153 #undef BLAS_DSCAL
154 #undef LAPACK_DPOTRF
155
156 #undef BLAS_ZTRSV
157 #undef BLAS_ZGEMV
158 #undef BLAS_ZTRSM
159 #undef BLAS_ZGEMM
160 #undef BLAS_ZHERK
161 #undef BLAS_ZGER
162 #undef BLAS_ZSCAL
163 #undef LAPACK_ZPOTRF
164
165
166
167 #define BLAS_DTRSV dtrsv_
168 #define BLAS_DGEMV dgemv_
169 #define BLAS_DTRSM dtrsm_
170 #define BLAS_DGEMM dgemm_
171 #define BLAS_DSYRK dsyrk_
172 #define BLAS_DGER dger_
173 #define BLAS_DSCAL dscal_
174 #define LAPACK_DPOTRF dpotrf_
175
176 #define BLAS_ZTRSV ztrsv_
177 #define BLAS_ZGEMV zgemv_
178 #define BLAS_ZTRSM ztrsm_
179 #define BLAS_ZGEMM zgemm_
180 #define BLAS_ZHERK zherk_
181 #define BLAS_ZGER zgeru_
182 #define BLAS_ZSCAL zscal_
183 #define LAPACK_ZPOTRF zpotrf_
184
185 #endif
/* ========================================================================== */ 145 186 /* ========================================================================== */
/* === BLAS and LAPACK integer arguments ==================================== */ 146 187 /* === BLAS and LAPACK integer arguments ==================================== */
/* ========================================================================== */ 147 188 /* ========================================================================== */
148 189
/* CHOLMOD can be compiled with -D'LONGBLAS=long' for the Sun Performance 149 190 /* CHOLMOD can be compiled with -D'LONGBLAS=long' for the Sun Performance
* Library, or -D'LONGBLAS=long long' for SGI's SCSL BLAS. This defines the 150 191 * Library, or -D'LONGBLAS=long long' for SGI's SCSL BLAS. This defines the
* integer used in the BLAS for the cholmod_l_* routines. 151 192 * integer used in the BLAS for the cholmod_l_* routines.
* 152 193 *
* The "int" version of CHOLMOD always uses the "int" version of the BLAS. 153 194 * The "int" version of CHOLMOD always uses the "int" version of the BLAS.
*/ 154 195 */
155 196
#if defined (LONGBLAS) && defined (LONG) 156 197 #if defined (LONGBLAS) && defined (LONG)
#define BLAS_INT LONGBLAS 157 198 #define BLAS_INT LONGBLAS
#else 158 199 #else
#define BLAS_INT int 159 200 #define BLAS_INT int
#endif 160 201 #endif
161 202
/* If the BLAS integer is smaller than the basic CHOLMOD integer, then we need 162 203 /* If the BLAS integer is smaller than the basic CHOLMOD integer, then we need
* to check for integer overflow when converting from one to the other. If 163 204 * to check for integer overflow when converting from one to the other. If
* any integer overflows, the externally-defined blas_ok variable is set to 164 205 * any integer overflows, the externally-defined blas_ok variable is set to
* FALSE. blas_ok should be set to TRUE before calling any BLAS_* macro. 165 206 * FALSE. blas_ok should be set to TRUE before calling any BLAS_* macro.
*/ 166 207 */
167 208
#define CHECK_BLAS_INT (sizeof (BLAS_INT) < sizeof (Int)) 168 209 #define CHECK_BLAS_INT (sizeof (BLAS_INT) < sizeof (Int))
#define EQ(K,k) (((BLAS_INT) K) == ((Int) k)) 169 210 #define EQ(K,k) (((BLAS_INT) K) == ((Int) k))
170 211
/* ========================================================================== */ 171 212 /* ========================================================================== */
/* === BLAS and LAPACK prototypes and macros ================================ */ 172 213 /* === BLAS and LAPACK prototypes and macros ================================ */
/* ========================================================================== */ 173 214 /* ========================================================================== */
174 215
void BLAS_DGEMV (char *trans, BLAS_INT *m, BLAS_INT *n, double *alpha, 175 216 void BLAS_DGEMV (char *trans, BLAS_INT *m, BLAS_INT *n, double *alpha,
double *A, BLAS_INT *lda, double *X, BLAS_INT *incx, double *beta, 176 217 double *A, BLAS_INT *lda, double *X, BLAS_INT *incx, double *beta,
double *Y, BLAS_INT *incy) ; 177 218 double *Y, BLAS_INT *incy) ;
178 219
#define BLAS_dgemv(trans,m,n,alpha,A,lda,X,incx,beta,Y,incy) \ 179 220 #define BLAS_dgemv(trans,m,n,alpha,A,lda,X,incx,beta,Y,incy) \
{ \ 180 221 { \
BLAS_INT M = m, N = n, LDA = lda, INCX = incx, INCY = incy ; \ 181 222 BLAS_INT M = m, N = n, LDA = lda, INCX = incx, INCY = incy ; \
if (CHECK_BLAS_INT) \ 182 223 if (CHECK_BLAS_INT) \
{ \ 183 224 { \
blas_ok &= EQ (M,m) && EQ (N,n) && EQ (LDA,lda) && EQ (INCX,incx) \ 184 225 blas_ok &= EQ (M,m) && EQ (N,n) && EQ (LDA,lda) && EQ (INCX,incx) \
&& EQ (INCY,incy) ; \ 185 226 && EQ (INCY,incy) ; \
} \ 186 227 } \
if (blas_ok) \ 187 228 if (blas_ok) \
{ \ 188 229 { \
BLAS_DGEMV (trans, &M, &N, alpha, A, &LDA, X, &INCX, beta, Y, &INCY) ; \ 189 230 BLAS_DGEMV (trans, &M, &N, alpha, A, &LDA, X, &INCX, beta, Y, &INCY) ; \
} \ 190 231 } \
} 191 232 }
192 233
void BLAS_ZGEMV (char *trans, BLAS_INT *m, BLAS_INT *n, double *alpha, 193 234 void BLAS_ZGEMV (char *trans, BLAS_INT *m, BLAS_INT *n, double *alpha,
double *A, BLAS_INT *lda, double *X, BLAS_INT *incx, double *beta, 194 235 double *A, BLAS_INT *lda, double *X, BLAS_INT *incx, double *beta,
double *Y, BLAS_INT *incy) ; 195 236 double *Y, BLAS_INT *incy) ;
196 237
#define BLAS_zgemv(trans,m,n,alpha,A,lda,X,incx,beta,Y,incy) \ 197 238 #define BLAS_zgemv(trans,m,n,alpha,A,lda,X,incx,beta,Y,incy) \
{ \ 198 239 { \
BLAS_INT M = m, N = n, LDA = lda, INCX = incx, INCY = incy ; \ 199 240 BLAS_INT M = m, N = n, LDA = lda, INCX = incx, INCY = incy ; \
if (CHECK_BLAS_INT) \ 200 241 if (CHECK_BLAS_INT) \
{ \ 201 242 { \
blas_ok &= EQ (M,m) && EQ (N,n) && EQ (LDA,lda) && EQ (INCX,incx) \ 202 243 blas_ok &= EQ (M,m) && EQ (N,n) && EQ (LDA,lda) && EQ (INCX,incx) \
&& EQ (INCY,incy) ; \ 203 244 && EQ (INCY,incy) ; \
} \ 204 245 } \
if (blas_ok) \ 205 246 if (blas_ok) \
{ \ 206 247 { \
BLAS_ZGEMV (trans, &M, &N, alpha, A, &LDA, X, &INCX, beta, Y, &INCY) ; \ 207 248 BLAS_ZGEMV (trans, &M, &N, alpha, A, &LDA, X, &INCX, beta, Y, &INCY) ; \
} \ 208 249 } \
} 209 250 }
210 251
void BLAS_DTRSV (char *uplo, char *trans, char *diag, BLAS_INT *n, double *A, 211 252 void BLAS_DTRSV (char *uplo, char *trans, char *diag, BLAS_INT *n, double *A,
BLAS_INT *lda, double *X, BLAS_INT *incx) ; 212 253 BLAS_INT *lda, double *X, BLAS_INT *incx) ;
213 254
#define BLAS_dtrsv(uplo,trans,diag,n,A,lda,X,incx) \ 214 255 #define BLAS_dtrsv(uplo,trans,diag,n,A,lda,X,incx) \
{ \ 215 256 { \
BLAS_INT N = n, LDA = lda, INCX = incx ; \ 216 257 BLAS_INT N = n, LDA = lda, INCX = incx ; \
if (CHECK_BLAS_INT) \ 217 258 if (CHECK_BLAS_INT) \
{ \ 218 259 { \
blas_ok &= EQ (N,n) && EQ (LDA,lda) && EQ (INCX,incx) ; \ 219 260 blas_ok &= EQ (N,n) && EQ (LDA,lda) && EQ (INCX,incx) ; \
} \ 220 261 } \
if (blas_ok) \ 221 262 if (blas_ok) \
{ \ 222 263 { \
BLAS_DTRSV (uplo, trans, diag, &N, A, &LDA, X, &INCX) ; \ 223 264 BLAS_DTRSV (uplo, trans, diag, &N, A, &LDA, X, &INCX) ; \
} \ 224 265 } \
} 225 266 }
226 267
void BLAS_ZTRSV (char *uplo, char *trans, char *diag, BLAS_INT *n, double *A, 227 268 void BLAS_ZTRSV (char *uplo, char *trans, char *diag, BLAS_INT *n, double *A,
BLAS_INT *lda, double *X, BLAS_INT *incx) ; 228 269 BLAS_INT *lda, double *X, BLAS_INT *incx) ;
229 270
#define BLAS_ztrsv(uplo,trans,diag,n,A,lda,X,incx) \ 230 271 #define BLAS_ztrsv(uplo,trans,diag,n,A,lda,X,incx) \
{ \ 231 272 { \
BLAS_INT N = n, LDA = lda, INCX = incx ; \ 232 273 BLAS_INT N = n, LDA = lda, INCX = incx ; \
if (CHECK_BLAS_INT) \ 233 274 if (CHECK_BLAS_INT) \
{ \ 234 275 { \
blas_ok &= EQ (N,n) && EQ (LDA,lda) && EQ (INCX,incx) ; \ 235 276 blas_ok &= EQ (N,n) && EQ (LDA,lda) && EQ (INCX,incx) ; \
} \ 236 277 } \
if (blas_ok) \ 237 278 if (blas_ok) \
{ \ 238 279 { \
BLAS_ZTRSV (uplo, trans, diag, &N, A, &LDA, X, &INCX) ; \ 239 280 BLAS_ZTRSV (uplo, trans, diag, &N, A, &LDA, X, &INCX) ; \
} \ 240 281 } \
} 241 282 }
242 283
void BLAS_DTRSM (char *side, char *uplo, char *transa, char *diag, BLAS_INT *m, 243 284 void BLAS_DTRSM (char *side, char *uplo, char *transa, char *diag, BLAS_INT *m,
BLAS_INT *n, double *alpha, double *A, BLAS_INT *lda, double *B, 244 285 BLAS_INT *n, double *alpha, double *A, BLAS_INT *lda, double *B,
BLAS_INT *ldb) ; 245 286 BLAS_INT *ldb) ;
246 287
#define BLAS_dtrsm(side,uplo,transa,diag,m,n,alpha,A,lda,B,ldb) \ 247 288 #define BLAS_dtrsm(side,uplo,transa,diag,m,n,alpha,A,lda,B,ldb) \
{ \ 248 289 { \
BLAS_INT M = m, N = n, LDA = lda, LDB = ldb ; \ 249 290 BLAS_INT M = m, N = n, LDA = lda, LDB = ldb ; \
if (CHECK_BLAS_INT) \ 250 291 if (CHECK_BLAS_INT) \
{ \ 251 292 { \
blas_ok &= EQ (M,m) && EQ (N,n) && EQ (LDA,lda) && EQ (LDB,ldb) ; \ 252 293 blas_ok &= EQ (M,m) && EQ (N,n) && EQ (LDA,lda) && EQ (LDB,ldb) ; \
} \ 253 294 } \
if (blas_ok) \ 254 295 if (blas_ok) \
{ \ 255 296 { \
BLAS_DTRSM (side, uplo, transa, diag, &M, &N, alpha, A, &LDA, B, &LDB);\ 256 297 BLAS_DTRSM (side, uplo, transa, diag, &M, &N, alpha, A, &LDA, B, &LDB);\
} \ 257 298 } \
} 258 299 }
259 300
void BLAS_ZTRSM (char *side, char *uplo, char *transa, char *diag, BLAS_INT *m, 260 301 void BLAS_ZTRSM (char *side, char *uplo, char *transa, char *diag, BLAS_INT *m,
BLAS_INT *n, double *alpha, double *A, BLAS_INT *lda, double *B, 261 302 BLAS_INT *n, double *alpha, double *A, BLAS_INT *lda, double *B,
BLAS_INT *ldb) ; 262 303 BLAS_INT *ldb) ;
263 304
#define BLAS_ztrsm(side,uplo,transa,diag,m,n,alpha,A,lda,B,ldb) \ 264 305 #define BLAS_ztrsm(side,uplo,transa,diag,m,n,alpha,A,lda,B,ldb) \
{ \ 265 306 { \
BLAS_INT M = m, N = n, LDA = lda, LDB = ldb ; \ 266 307 BLAS_INT M = m, N = n, LDA = lda, LDB = ldb ; \
if (CHECK_BLAS_INT) \ 267 308 if (CHECK_BLAS_INT) \
{ \ 268 309 { \
blas_ok &= EQ (M,m) && EQ (N,n) && EQ (LDA,lda) && EQ (LDB,ldb) ; \ 269 310 blas_ok &= EQ (M,m) && EQ (N,n) && EQ (LDA,lda) && EQ (LDB,ldb) ; \
} \ 270 311 } \
if (blas_ok) \ 271 312 if (blas_ok) \
{ \ 272 313 { \
BLAS_ZTRSM (side, uplo, transa, diag, &M, &N, alpha, A, &LDA, B, &LDB);\ 273 314 BLAS_ZTRSM (side, uplo, transa, diag, &M, &N, alpha, A, &LDA, B, &LDB);\
} \ 274 315 } \
} 275 316 }
276 317
void BLAS_DGEMM (char *transa, char *transb, BLAS_INT *m, BLAS_INT *n, 277 318 void BLAS_DGEMM (char *transa, char *transb, BLAS_INT *m, BLAS_INT *n,
BLAS_INT *k, double *alpha, double *A, BLAS_INT *lda, double *B, 278 319 BLAS_INT *k, double *alpha, double *A, BLAS_INT *lda, double *B,
BLAS_INT *ldb, double *beta, double *C, BLAS_INT *ldc) ; 279 320 BLAS_INT *ldb, double *beta, double *C, BLAS_INT *ldc) ;
280 321
#define BLAS_dgemm(transa,transb,m,n,k,alpha,A,lda,B,ldb,beta,C,ldc) \ 281 322 #define BLAS_dgemm(transa,transb,m,n,k,alpha,A,lda,B,ldb,beta,C,ldc) \
{ \ 282 323 { \
BLAS_INT M = m, N = n, K = k, LDA = lda, LDB = ldb, LDC = ldc ; \ 283 324 BLAS_INT M = m, N = n, K = k, LDA = lda, LDB = ldb, LDC = ldc ; \
if (CHECK_BLAS_INT) \ 284 325 if (CHECK_BLAS_INT) \
{ \ 285 326 { \
blas_ok &= EQ (M,m) && EQ (N,n) && EQ (K,k) && EQ (LDA,lda) \ 286 327 blas_ok &= EQ (M,m) && EQ (N,n) && EQ (K,k) && EQ (LDA,lda) \
&& EQ (LDB,ldb) && EQ (LDC,ldc) ; \ 287 328 && EQ (LDB,ldb) && EQ (LDC,ldc) ; \
} \ 288 329 } \
if (blas_ok) \ 289 330 if (blas_ok) \
{ \ 290 331 { \
BLAS_DGEMM (transa, transb, &M, &N, &K, alpha, A, &LDA, B, &LDB, beta, \ 291 332 BLAS_DGEMM (transa, transb, &M, &N, &K, alpha, A, &LDA, B, &LDB, beta, \
C, &LDC) ; \ 292 333 C, &LDC) ; \
} \ 293 334 } \
} 294 335 }
295 336
void BLAS_ZGEMM (char *transa, char *transb, BLAS_INT *m, BLAS_INT *n, 296 337 void BLAS_ZGEMM (char *transa, char *transb, BLAS_INT *m, BLAS_INT *n,
BLAS_INT *k, double *alpha, double *A, BLAS_INT *lda, double *B, 297 338 BLAS_INT *k, double *alpha, double *A, BLAS_INT *lda, double *B,
BLAS_INT *ldb, double *beta, double *C, BLAS_INT *ldc) ; 298 339 BLAS_INT *ldb, double *beta, double *C, BLAS_INT *ldc) ;
299 340
#define BLAS_zgemm(transa,transb,m,n,k,alpha,A,lda,B,ldb,beta,C,ldc) \ 300 341 #define BLAS_zgemm(transa,transb,m,n,k,alpha,A,lda,B,ldb,beta,C,ldc) \
{ \ 301 342 { \
BLAS_INT M = m, N = n, K = k, LDA = lda, LDB = ldb, LDC = ldc ; \ 302 343 BLAS_INT M = m, N = n, K = k, LDA = lda, LDB = ldb, LDC = ldc ; \
if (CHECK_BLAS_INT) \ 303 344 if (CHECK_BLAS_INT) \
{ \ 304 345 { \
blas_ok &= EQ (M,m) && EQ (N,n) && EQ (K,k) && EQ (LDA,lda) \ 305 346 blas_ok &= EQ (M,m) && EQ (N,n) && EQ (K,k) && EQ (LDA,lda) \
&& EQ (LDB,ldb) && EQ (LDC,ldc) ; \ 306 347 && EQ (LDB,ldb) && EQ (LDC,ldc) ; \
} \ 307 348 } \
if (blas_ok) \ 308 349 if (blas_ok) \
{ \ 309 350 { \
BLAS_ZGEMM (transa, transb, &M, &N, &K, alpha, A, &LDA, B, &LDB, beta, \ 310 351 BLAS_ZGEMM (transa, transb, &M, &N, &K, alpha, A, &LDA, B, &LDB, beta, \
C, &LDC) ; \ 311 352 C, &LDC) ; \
} \ 312 353 } \
} 313 354 }
314 355
void BLAS_DSYRK (char *uplo, char *trans, BLAS_INT *n, BLAS_INT *k, 315 356 void BLAS_DSYRK (char *uplo, char *trans, BLAS_INT *n, BLAS_INT *k,
double *alpha, double *A, BLAS_INT *lda, double *beta, double *C, 316 357 double *alpha, double *A, BLAS_INT *lda, double *beta, double *C,
BLAS_INT *ldc) ; 317 358 BLAS_INT *ldc) ;
318 359
#define BLAS_dsyrk(uplo,trans,n,k,alpha,A,lda,beta,C,ldc) \ 319 360 #define BLAS_dsyrk(uplo,trans,n,k,alpha,A,lda,beta,C,ldc) \
{ \ 320 361 { \
BLAS_INT N = n, K = k, LDA = lda, LDC = ldc ; \ 321 362 BLAS_INT N = n, K = k, LDA = lda, LDC = ldc ; \
if (CHECK_BLAS_INT) \ 322 363 if (CHECK_BLAS_INT) \
{ \ 323 364 { \
blas_ok &= EQ (N,n) && EQ (K,k) && EQ (LDA,lda) && EQ (LDC,ldc) ; \ 324 365 blas_ok &= EQ (N,n) && EQ (K,k) && EQ (LDA,lda) && EQ (LDC,ldc) ; \
} \ 325 366 } \
if (blas_ok) \ 326 367 if (blas_ok) \
{ \ 327 368 { \
BLAS_DSYRK (uplo, trans, &N, &K, alpha, A, &LDA, beta, C, &LDC) ; \ 328 369 BLAS_DSYRK (uplo, trans, &N, &K, alpha, A, &LDA, beta, C, &LDC) ; \
} \ 329 370 } \
} \ 330 371 } \
331 372
void BLAS_ZHERK (char *uplo, char *trans, BLAS_INT *n, BLAS_INT *k, 332 373 void BLAS_ZHERK (char *uplo, char *trans, BLAS_INT *n, BLAS_INT *k,
double *alpha, double *A, BLAS_INT *lda, double *beta, double *C, 333 374 double *alpha, double *A, BLAS_INT *lda, double *beta, double *C,
BLAS_INT *ldc) ; 334 375 BLAS_INT *ldc) ;
335 376
#define BLAS_zherk(uplo,trans,n,k,alpha,A,lda,beta,C,ldc) \ 336 377 #define BLAS_zherk(uplo,trans,n,k,alpha,A,lda,beta,C,ldc) \
{ \ 337 378 { \
BLAS_INT N = n, K = k, LDA = lda, LDC = ldc ; \ 338 379 BLAS_INT N = n, K = k, LDA = lda, LDC = ldc ; \
if (CHECK_BLAS_INT) \ 339 380 if (CHECK_BLAS_INT) \
{ \ 340 381 { \
blas_ok &= EQ (N,n) && EQ (K,k) && EQ (LDA,lda) && EQ (LDC,ldc) ; \ 341 382 blas_ok &= EQ (N,n) && EQ (K,k) && EQ (LDA,lda) && EQ (LDC,ldc) ; \
} \ 342 383 } \
if (blas_ok) \ 343 384 if (blas_ok) \
{ \ 344 385 { \
BLAS_ZHERK (uplo, trans, &N, &K, alpha, A, &LDA, beta, C, &LDC) ; \ 345 386 BLAS_ZHERK (uplo, trans, &N, &K, alpha, A, &LDA, beta, C, &LDC) ; \
} \ 346 387 } \
} \ 347 388 } \
348 389
void LAPACK_DPOTRF (char *uplo, BLAS_INT *n, double *A, BLAS_INT *lda, 349 390 void LAPACK_DPOTRF (char *uplo, BLAS_INT *n, double *A, BLAS_INT *lda,
BLAS_INT *info) ; 350 391 BLAS_INT *info) ;
351 392
#define LAPACK_dpotrf(uplo,n,A,lda,info) \ 352 393 #define LAPACK_dpotrf(uplo,n,A,lda,info) \
{ \ 353 394 { \
BLAS_INT N = n, LDA = lda, INFO = 1 ; \ 354 395 BLAS_INT N = n, LDA = lda, INFO = 1 ; \
if (CHECK_BLAS_INT) \ 355 396 if (CHECK_BLAS_INT) \
{ \ 356 397 { \
blas_ok &= EQ (N,n) && EQ (LDA,lda) ; \ 357 398 blas_ok &= EQ (N,n) && EQ (LDA,lda) ; \
} \ 358 399 } \
if (blas_ok) \ 359 400 if (blas_ok) \
{ \ 360 401 { \
LAPACK_DPOTRF (uplo, &N, A, &LDA, &INFO) ; \ 361 402 LAPACK_DPOTRF (uplo, &N, A, &LDA, &INFO) ; \
} \ 362 403 } \
info = INFO ; \ 363 404 info = INFO ; \
} 364 405 }
365 406
void LAPACK_ZPOTRF (char *uplo, BLAS_INT *n, double *A, BLAS_INT *lda, 366 407 void LAPACK_ZPOTRF (char *uplo, BLAS_INT *n, double *A, BLAS_INT *lda,
BLAS_INT *info) ; 367 408 BLAS_INT *info) ;
368 409
#define LAPACK_zpotrf(uplo,n,A,lda,info) \ 369 410 #define LAPACK_zpotrf(uplo,n,A,lda,info) \
{ \ 370 411 { \
BLAS_INT N = n, LDA = lda, INFO = 1 ; \ 371 412 BLAS_INT N = n, LDA = lda, INFO = 1 ; \
if (CHECK_BLAS_INT) \ 372 413 if (CHECK_BLAS_INT) \
{ \ 373 414 { \
blas_ok &= EQ (N,n) && EQ (LDA,lda) ; \ 374 415 blas_ok &= EQ (N,n) && EQ (LDA,lda) ; \
} \ 375 416 } \
if (blas_ok) \ 376 417 if (blas_ok) \
{ \ 377 418 { \
LAPACK_ZPOTRF (uplo, &N, A, &LDA, &INFO) ; \ 378 419 LAPACK_ZPOTRF (uplo, &N, A, &LDA, &INFO) ; \
} \ 379 420 } \
info = INFO ; \ 380 421 info = INFO ; \
} 381 422 }
382 423
/* ========================================================================== */ 383 424 /* ========================================================================== */
384 425
1 1
include $(BTREE)/Make.inc 2 2 include $(BTREE)/Make.inc
3 3
programs = test_fac test_matinv test_specfunc \ 4 4 programs = test_fac$(exext) test_matinv$(exext) test_specfunc$(exext) \
test_det test_matcon test_matev test_sparse test_inter1d \ 5 5 test_det$(exext) test_matcon$(exext) test_matev$(exext) test_sparse$(exext) test_inter1d$(exext) \
test_inter2d test_inter3d test_akima test_lsp test_muller \ 6 6 test_inter2d$(exext) test_inter3d$(exext) test_akima$(exext) test_lsp$(exext) test_muller$(exext) \
test_integ test_bsyn test_bsjn test_bskn test_bsin 7 7 test_integ$(exext) test_bsyn$(exext) test_bsjn$(exext) test_bskn$(exext) test_bsin$(exext)
8 8
prog:$(programs) 9 9 prog:$(programs)
10 10
clean: 11 11 clean:
rm -f {*.o,*.oo,*.ipo,*.a,*.mod} 12 12 rm -f {*.o,*.oo,*.ipo,*.a,*.mod}
rm -f $(programs) 13 13 rm -f $(programs)
14 14
%: %.o 15 15 %$(exext): %.o
$(LINK) $(LINKFLAGS) $< init_random_seed.o $(LINKFVN) -o $@ 16 16 $(LINK) $(LINKFLAGS) $< init_random_seed.o $(LINKFVN) -o $@
17 17
%.o: %.f90 18 18 %.o: %.f90
$(F95) $(F95FLAGS) -c $< 19 19 $(F95) $(F95FLAGS) -c $< -o $@
20 20
fvn_test/test_akima.f90
program akima 1 1 program akima
use fvn 2 2 use fvn
implicit none 3 3 implicit none
integer :: nbpoints,nppoints,i 4 4 integer :: nbpoints,nppoints,i
real(8),dimension(:),allocatable :: x_d,y_d,breakpoints_d 5 5 real(8),dimension(:),allocatable :: x_d,y_d,breakpoints_d
real(8),dimension(:,:),allocatable :: coeff_fvn_d 6 6 real(8),dimension(:,:),allocatable :: coeff_fvn_d
real(8) :: xstep_d,xp_d,ty_d,fvn_y_d 7 7 real(8) :: xstep_d,xp_d,ty_d,fvn_y_d
open(2,file='fvn_akima_double.dat') 8 8 open(2,file='fvn_akima_double.dat')
open(3,file='fvn_akima_breakpoints_double.dat') 9 9 open(3,file='fvn_akima_breakpoints_double.dat')
nbpoints=30 10 10 nbpoints=30
allocate(x_d(nbpoints)) 11 11 allocate(x_d(nbpoints))
allocate(y_d(nbpoints)) 12 12 allocate(y_d(nbpoints))
allocate(breakpoints_d(nbpoints)) 13 13 allocate(breakpoints_d(nbpoints))
allocate(coeff_fvn_d(4,nbpoints)) 14 14 allocate(coeff_fvn_d(4,nbpoints))
xstep_d=20./dfloat(nbpoints) 15 15 xstep_d=20./dfloat(nbpoints)
do i=1,nbpoints 16 16 do i=1,nbpoints
x_d(i)=-10.+dfloat(i)*xstep_d 17 17 x_d(i)=-10.+dfloat(i)*xstep_d
y_d(i)=dsin(x_d(i)) 18 18 y_d(i)=dsin(x_d(i))
write(3,44) x_d(i),y_d(i) 19 19 write(3,44) x_d(i),y_d(i)
end do 20 20 end do
close(3) 21 21 close(3)
call fvn_akima(nbpoints,x_d,y_d,breakpoints_d,coeff_fvn_d) 22 22 call fvn_akima(nbpoints,x_d,y_d,breakpoints_d,coeff_fvn_d)
nppoints=1000 23 23 nppoints=1000
xstep_d=22./dfloat(nppoints) 24 24 xstep_d=22./dfloat(nppoints)
do i=1,nppoints 25 25 do i=1,nppoints
xp_d=-11.+dfloat(i)*xstep_d 26 26 xp_d=-11.+dfloat(i)*xstep_d
ty_d=dsin(xp_d) 27 27 ty_d=dsin(xp_d)
fvn_y_d=fvn_spline_eval(xp_d,nbpoints-1,breakpoints_d,coeff_fvn_d) 28 28 fvn_y_d=fvn_spline_eval(xp_d,nbpoints-1,breakpoints_d,coeff_fvn_d)
write(2,44) xp_d,ty_d,fvn_y_d 29 29 write(2,44) xp_d,ty_d,fvn_y_d
end do 30 30 end do
close(2) 31 31 close(2)
deallocate(coeff_fvn_d,breakpoints_d,y_d,x_d) 32 32 deallocate(coeff_fvn_d,breakpoints_d,y_d,x_d)
write(*,*) "All done, plot results with gnuplot using command :" 33 33 write(6,*) "All done, plot results with gnuplot using command :"
write(*,*) "pl 'fvn_akima_double.dat' u 1:2 w l,'fvn_akima_breakpoints_double.dat' w p" 34 34 write(6,*) "pl 'fvn_akima_double.dat' u 1:2 w l,'fvn_akima_breakpoints_double.dat' w p"
35 35
44 FORMAT(4(1X,1PE22.14)) 36 36 44 FORMAT(4(1X,1PE22.14))
end program 37 37 end program
38 38
fvn_test/test_sparse.f90
program test_sparse 1 1 program test_sparse
use fvn 2 2 use fvn
implicit none 3 3 implicit none
integer(4), parameter :: nz=12 4 4 integer(4), parameter :: nz=12
integer(4), parameter :: n=5 5 5 integer(4), parameter :: n=5
real(8),dimension(nz) :: A 6 6 real(8),dimension(nz) :: A
real(8),dimension(n,n) :: As 7 7 real(8),dimension(n,n) :: As
integer(4),dimension(nz) :: Ti,Tj 8 8 integer(4),dimension(nz) :: Ti,Tj
real(8),dimension(n) :: B,x 9 9 real(8),dimension(n) :: B,x
integer(4) :: status,i 10 10 integer(4) :: status,i
! Description of the matrix in triplet form 11 11 ! Description of the matrix in triplet form
A = (/ 2.,3.,3.,-1.,4.,4.,-3.,1.,2.,2.,6.,1. /) 12 12 A = (/ 2.,3.,3.,-1.,4.,4.,-3.,1.,2.,2.,6.,1. /)
B = (/ 8., 45., -3., 3., 19./) 13 13 B = (/ 8., 45., -3., 3., 19./)
Ti = (/ 1,2,1,3,5,2,3,4,5,3,2,5 /) 14 14 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 /) 15 15 Tj = (/ 1,1,2,2,2,3,3,3,3,4,5,5 /)
16 16
! Reconstruction of the matrix in standard form 17 17 ! Reconstruction of the matrix in standard form
! just needed for printing the matrix here 18 18 ! just needed for printing the matrix here
As=0. 19 19 As=0.
do i=1,nz 20 20 do i=1,nz
As(Ti(i),Tj(i))=A(i) 21 21 As(Ti(i),Tj(i))=A(i)
end do 22 22 end do
write(*,*) "Matrix in standard representation :" 23 23 write(*,*) "Matrix in standard representation :"
do i=1,5 24 24 do i=1,5
write(*,'(5f8.4)') As(i,:) 25 25 write(*,'(5f8.4)') As(i,:)
end do 26 26 end do
write(*,*) 27 27 write(*,*)
write(*,'("Right hand side :",5f8.4)') B 28 28 write(*,'("Right hand side :",5f8.4)') B
29 29
!specific routine that will be used here 30 30 !specific routine that will be used here
!call fvn_di_sparse_solve(n,nz,A,Ti,Tj,B,x,status) 31 31 !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) 32 32 call fvn_di_sparse_solve(n,nz,A,Ti,Tj,B,x,status)
write(*,'("Solution :",5f8.4)') x 33 33 write(*,'("Solution :",5f8.4)') x
write(*,'("Product matrix Solution :",5f8.4)') matmul(As,x) 34 34 write(*,'("Product matrix Solution :",5f8.4)') matmul(As,x)