diff --git a/ramsurflib.c b/ramsurflib.c index e03df71..5406c04 100644 --- a/ramsurflib.c +++ b/ramsurflib.c @@ -370,7 +370,7 @@ void deriv(int n, float sig, double alp, dh1[0][1]=sig*0.5; double exp1 = -0.5; dh2[0][0]=alp; - dh2[0][1]=alp; + dh2[0][1]=0.; double exp2=-1.0; dh3[0][0]=-2.0*nu; dh3[0][1]=0.; @@ -830,7 +830,7 @@ void matrc(size_t mz, size_t mp, int const nz, int const np, int const iz, float static void updat( ramsurf_t const* rsurf, size_t *profl_index, size_t mr, size_t mz, size_t mp, int nz, int np, int *iz, int *ib, float dr, float dz, float omega, float k0, - float r, float *rp, float rs, + float r, float *rp, float *rs, float rb[mr], float zb[mr] , float cw[mz], float cb[mz], float rhob[mz], float attn[mz], float alpw[mz], float alpb[mz], fcomplex ksq[mz], fcomplex ksqw[mz], fcomplex ksqb[mz], @@ -871,8 +871,9 @@ void updat( ramsurf_t const* rsurf, size_t *profl_index, size_t mr, size_t mz, s // // Turn off the stability constraints. // - if(r>=rs) { + if(r>=*rs) { int ns=0; + *rs=2.f*(rsurf->rmax); epade(mp, np, ns, 1, k0, dr, (float (*)[2])pd1, (float (*)[2])pd2); matrc(mz, mp, nz, np, *iz, dz, k0, rhob, alpw, alpb, (float (*)[2])ksq, (float (*)[2])ksqw, (float (*)[2])ksqb, f1, f2, f3, r1, r2, r3, s1, s2, s3, (float (*)[2])pd1, (float (*)[2])pd2, *izsrf); @@ -1118,6 +1119,13 @@ int ramsurf(ramsurf_t const* rsurf, int * lz, float *** ogrid, FILE *fdline) fix_zmax(&lzmax, rsurf->rhob); fix_zmax(&lzmax, rsurf->attn); +#ifdef __SSE3__ + // this sets the "flush to zero" and "denormals are zero" bits (15,6) in the MXCSR register. + // otherwise we get very bad timings... + int flush_mode = _MM_GET_FLUSH_ZERO_MODE( ); + _MM_SET_FLUSH_ZERO_MODE(_MM_FLUSH_ZERO_ON); +#endif + size_t mr = rsurf->mr, mz = lzmax / rsurf->dz + 2.5f, mp = rsurf->np; @@ -1177,7 +1185,7 @@ int ramsurf(ramsurf_t const* rsurf, int * lz, float *** ogrid, FILE *fdline) // while (r < rsurf->rmax) { updat(rsurf, &profl_index, mr, mz, mp, nz, np, &iz, &ib, dr, dz, omega, k0, r, - &rp, rs, *rb, *zb, *cw, *cb, *rhob, *attn, *alpw, *alpb, *ksq, *ksqw, *ksqb, *f1, *f2, *f3, + &rp, &rs, *rb, *zb, *cw, *cb, *rhob, *attn, *alpw, *alpb, *ksq, *ksqw, *ksqb, *f1, *f2, *f3, *r1, *r2, *r3, *s1, *s2, *s3, *pd1, *pd2, *rsrf, *zsrf, &izsrf, &isrf); solve(mz, mp, nz, np, *u, *r1, *r3, *s1, *s2, *s3); r=r+dr; @@ -1191,6 +1199,9 @@ int ramsurf(ramsurf_t const* rsurf, int * lz, float *** ogrid, FILE *fdline) // deallocation step free(oscratch); +#ifdef __SSE3__ + _MM_SET_FLUSH_ZERO_MODE(flush_mode); +#endif return errorCode; } diff --git a/tests/Makefile.am b/tests/Makefile.am index 6ccb019..93dbdb2 100644 --- a/tests/Makefile.am +++ b/tests/Makefile.am @@ -1,6 +1,7 @@ -check_PROGRAMS=ramcmp +check_PROGRAMS=ramcmp ramdump ramcmp_SOURCES=ramcmp.f +ramdump_SOURCES=ramdump.f TESTS=decreasing_slop.sh\ flat.sh\ @@ -9,7 +10,8 @@ TESTS=decreasing_slop.sh\ flat_pade.sh\ flat_limit.sh\ increasing_slop.sh\ - flat_2ssps.sh + flat_2ssps.sh\ + intensive.sh EXTRA_DIST=decreasing_slop.test\ flat.test\ @@ -18,18 +20,20 @@ EXTRA_DIST=decreasing_slop.test\ flat_pade.test\ flat_limit.test\ increasing_slop.test\ - flat_2ssps.test + flat_2ssps.test\ + intensive.test SUFFIXES=.test .sh .test.sh:$(srcdir)/Makefile.am echo '#!/bin/sh' > $@ echo 'set -e' >> $@ - echo 'chmod -R u+rwx `basename $<`.d || true' >> $@ - echo 'rm -rf `basename $<`.d' >> $@ - echo 'cp -r $(srcdir)/`basename $<` `basename $<`.d' >> $@ - echo 'chmod -R u+rwx `basename $<`.d' >> $@ - echo 'cd `basename $<`.d' >> $@ + echo 'TARGET="`basename $<`"' >> $@ + echo 'test ! -d $$TARGET.d || chmod -R u+rwx $$TARGET.d' >> $@ + echo 'rm -rf $$TARGET' >> $@ + echo 'cp -r $(srcdir)/$$TARGET $$TARGET.d' >> $@ + echo 'chmod -R u+rwx $$TARGET.d' >> $@ + echo 'cd $$TARGET.d' >> $@ echo '$(abs_top_builddir)/ramsurf1.5' >> $@ echo 'mv tl.grid ref.grid' >> $@ echo 'mv tl.line ref.line' >> $@ diff --git a/tests/intensive.test/ram.in b/tests/intensive.test/ram.in new file mode 100644 index 0000000..1368e48 --- /dev/null +++ b/tests/intensive.test/ram.in @@ -0,0 +1,24 @@ +header line +2048.00 5 1000 1.0 frequency (Hz), source depth (m), tl.line depth (m) +400000. 10 40 domain range (m), computation step (m), mesure step (comp.step) +1000 0.2 10 1000 domain depth (m), computation step (m), mesure step (comp.step), last mesure depth (m) +1500.000000 4 1 0.000000 reference celerity, #PADE, #stability constraints, radius of constraint +0.0 0.0 surface depth (range, depth) (m) +-1 -1 + 0.0000 1000. bathymetry +400000.000 900 +-1 -1 + 0.0000 1500 + 333.000 1500 + 334.000 1510 + 666.000 1510 + 667.000 1500 + 1000.000 1500 +-1 -1 +0.0000000 3000.0000000 +-1 -1 +0.0000000 2. +-1 -1 +0.0000000 0.900 +21 100.00 +-1 -1 diff --git a/tests/ramcmp.f b/tests/ramcmp.f index 8277f26..06dd1c6 100644 --- a/tests/ramcmp.f +++ b/tests/ramcmp.f @@ -1,7 +1,7 @@ program ramdec real tl(1000,1000) real rtl(1000,1000) - threshold = 6.8 + threshold = .5 open(unit=2,status='unknown',file='ref.grid', form='unformatted') open(unit=3,status='unknown',file='tl.grid', form='unformatted') read(2)nry diff --git a/tests/ramdump.f b/tests/ramdump.f new file mode 100644 index 0000000..48f9e4f --- /dev/null +++ b/tests/ramdump.f @@ -0,0 +1,15 @@ + program ramdec + real tl(1000,1000) + CHARACTER(len=256) :: arg + threshold = .3 + eps = 1e-10 + CALL get_command_argument(1, arg) + open(unit=2,status='unknown',file=arg, form='unformatted') + read(2)nry + nx=1 + 5 read(2,end=6)(tl(j,nx),j=1,nry) + write(*,*)(tl(j,nx) ,j=1,nry) + nx=nx+1 + go to 5 + 6 end +