Skip to content
Snippets Groups Projects
Commit 70267546 authored by John Davis Tingle's avatar John Davis Tingle
Browse files

Delete _livermorec.c_

parent 38a89836
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
Livermore Tests/_livermorec.c_ deleted 100644 → 0
+ 0
− 1435
View file @ 38a89836
/* Livermore Loops coded in C Latest File Modification 20 Oct 92,
* by Tim Peters, Kendall Square Res. Corp. tim@ksr.com, ksr!tim@uunet.uu.net
* SUBROUTINE KERNEL( TK) replaces the Fortran routine in LFK Test program.
************************************************************************
* *
* KERNEL executes 24 samples of "C" computation *
* *
* TK(1) - total cpu time to execute only the 24 kernels.*
* TK(2) - total Flops executed by the 24 Kernels *
* *
************************************************************************
* *
* L. L. N. L. " C " K E R N E L S: M F L O P S *
* *
* These kernels measure " C " numerical computation *
* rates for a spectrum of cpu-limited computational *
* structures or benchmarks. Mathematical through-put *
* is measured in units of millions of floating-point *
* operations executed per second, called Megaflops/sec. *
* *
* Fonzi's Law: There is not now and there never will be a language *
* in which it is the least bit difficult to write *
* bad programs. *
* F.H.MCMAHON 1972 *
************************************************************************
*Originally from Greg Astfalk, AT&T, P.O.Box 900, Princeton, NJ. 08540*
* by way of Frank McMahon (LLNL). *
* *
* REFERENCE *
* *
* F.H.McMahon, The Livermore Fortran Kernels: *
* A Computer Test Of The Numerical Performance Range, *
* Lawrence Livermore National Laboratory, *
* Livermore, California, UCRL-53745, December 1986. *
* *
* from: National Technical Information Service *
* U.S. Department of Commerce *
* 5285 Port Royal Road *
* Springfield, VA. 22161 *
* *
* Changes made to correct many array subscripting problems, *
* make more readable (added #define's), include the original *
* FORTRAN versions of the runs as comments, and make more *
* portable by Kelly O'Hair (LLNL) and Chuck Rasbold (LLNL). *
* *
************************************************************************
*/
#include <stdio.h>
#include <math.h>
/* Type-specifiers for function declarations */
#ifdef CRAY
#define CALLED_FROM_FORTRAN fortran
#define FORTRAN_FUNCTION fortran
#else
#define CALLED_FROM_FORTRAN
#define FORTRAN_FUNCTION extern
#endif
/* These are the names to be used for the functions */
#ifdef CRAY
#define KERNEL kernel
#define TRACE trace
#define SPACE space
#define TEST test
#define TRACK track
#else
#define KERNEL kernel_
#define TRACE trace_
#define SPACE space_
#define TEST test_
#define TRACK track_
#endif
/* Type-specifiers for the structs that map to common */
#ifdef CRAY
#define COMMON_BLOCK fortran
#else
#define COMMON_BLOCK extern
#endif
/* Names of structs (or common blocks) */
#ifdef CRAY
#define ALPHA alpha
#define BETA beta
#define SPACES spaces
#define SPACER spacer
#define SPACE0 space0
#define SPACEI spacei
#define ISPACE ispace
#define SPACE1 space1
#define SPACE2 space2
#else
#define ALPHA alpha_
#define BETA beta_
#define SPACES spaces_
#define SPACER spacer_
#define SPACE0 space0_
#define SPACEI spacei_
#define ISPACE ispace_
#define SPACE1 space1_
#define SPACE2 space2_
#endif
/* Declare the four fortran functions called */
FORTRAN_FUNCTION TRACE();
FORTRAN_FUNCTION SPACE();
FORTRAN_FUNCTION TEST();
FORTRAN_FUNCTION TRACK();
/* Define the structs or COMMON BLOCKS */
COMMON_BLOCK struct {
long Mk;
long Ik;
long Im;
long Ml;
long Il;
long Mruns;
long Nruns;
long Jr;
long Npfs[47][3][8];
} ALPHA ;
#define mk ALPHA.Mk
#define ik ALPHA.Ik
#define im ALPHA.Im
#define ml ALPHA.Ml
#define il ALPHA.Il
#define mruns ALPHA.Mruns;
#define nruns ALPHA.Nruns;
#define jr ALPHA.Jr
#define npfs ALPHA.Npfs
COMMON_BLOCK struct {
double Tic;
double Times[47][3][8];
double See[3][8][3][5];
double Terrs[47][3][8];
double Csums[47][3][8];
double Fopn[47][3][8];
double Dos[47][3][8];
} BETA ;
#define tic BETA.Tic
#define times BETA.Times
#define see BETA.See
#define terrs BETA.Terrs
#define csums BETA.Csums
#define fopn BETA.Fopn
#define dos BETA.Dos
COMMON_BLOCK struct {
long Ion;
long J5;
long K2;
long K3;
long MULTI;
long Laps;
long Loop;
long M;
long Kr;
long It;
long N13h;
long Ibuf;
long Npass;
long Nfail;
long N;
long N1;
long N2;
long N13;
long N213;
long N813;
long N14;
long N16;
long N416;
long N21;
long Nt1;
long Nt2;
} SPACES ;
#define ion SPACES.Ion
#define j5 SPACES.J5
#define k2 SPACES.K2
#define k3 SPACES.K3
#define multi SPACES.MULTI
#define laps SPACES.Laps
#define loop SPACES.Loop
#define m SPACES.M
#define kr SPACES.Kr
#define it SPACES.It
#define n13h SPACES.N13h
#define ibuf SPACES.Ibuf
#define npass SPACES.Npass
#define nfail SPACES.Nfail
#define n SPACES.N
#define n1 SPACES.N1
#define n2 SPACES.N2
#define n13 SPACES.N13
#define n213 SPACES.N213
#define n813 SPACES.N813
#define n14 SPACES.N14
#define n16 SPACES.N16
#define n416 SPACES.N416
#define n21 SPACES.N21
#define nt1 SPACES.Nt1
#define nt2 SPACES.Nt2
COMMON_BLOCK struct {
double A11;
double A12;
double A13;
double A21;
double A22;
double A23;
double A31;
double A32;
double A33;
double Ar;
double Br;
double C0;
double Cr;
double Di;
double Dk;
double Dm22;
double Dm23;
double Dm24;
double Dm25;
double Dm26;
double Dm27;
double Dm28;
double Dn;
double E3;
double E6;
double Expmax;
double Flx;
double Q;
double Qa;
double R;
double Ri;
double S;
double Scale;
double Sig;
double Stb5;
double T;
double Xnc;
double Xnei;
double Xnm;
} SPACER ;
#define a11 SPACER.A11
#define a12 SPACER.A12
#define a13 SPACER.A13
#define a21 SPACER.A21
#define a22 SPACER.A22
#define a23 SPACER.A23
#define a31 SPACER.A31
#define a32 SPACER.A32
#define a33 SPACER.A33
#define ar SPACER.Ar
#define br SPACER.Br
#define c0 SPACER.C0
#define cr SPACER.Cr
#define di SPACER.Di
#define dk SPACER.Dk
#define dm22 SPACER.Dm22
#define dm23 SPACER.Dm23
#define dm24 SPACER.Dm24
#define dm25 SPACER.Dm25
#define dm26 SPACER.Dm26
#define dm27 SPACER.Dm27
#define dm28 SPACER.Dm28
#define dn SPACER.Dn
#define e3 SPACER.E3
#define e6 SPACER.E6
#define expmax SPACER.Expmax
#define flx SPACER.Flx
#define q SPACER.Q
#define qa SPACER.Qa
#define r SPACER.R
#define ri SPACER.Ri
#define s SPACER.S
#define scale SPACER.Scale
#define sig SPACER.Sig
#define stb5 SPACER.Stb5
#define t SPACER.T
#define xnc SPACER.Xnc
#define xnei SPACER.Xnei
#define xnm SPACER.Xnm
COMMON_BLOCK struct {
double Time[47];
double Csum[47];
double Ww[47];
double Wt[47];
double Ticks;
double Fr[9];
double Terr1[47];
double Sumw[7];
double Start;
double Skale[47];
double Bias[47];
double Ws[95];
double Total[47];
double Flopn[47];
long Iq[7];
long Npf;
long Npfs1[47];
} SPACE0 ;
#define time SPACE0.Time
#define csum SPACE0.Csum
#define ww SPACE0.Ww
#define wt SPACE0.Wt
#define ticks SPACE0.Ticks
#define fr SPACE0.Fr
#define terr1 SPACE0.Terr1
#define sumw SPACE0.Sumw
#define start SPACE0.Start
#define skale SPACE0.Skale
#define bias SPACE0.Bias
#define ws SPACE0.Ws
#define total SPACE0.Total
#define flopn SPACE0.Flopn
#define iq SPACE0.Iq
#define npf SPACE0.Npf
#define npfs1 SPACE0.Npfs1
COMMON_BLOCK struct {
double Wtp[3];
long Mult[3];
long Ispan[3][47];
long Ipass[3][47];
} SPACEI ;
#define wtp SPACEI.Wtp
#define mult SPACEI.Mult
#define ispan SPACEI.Ispan
#define ipass SPACEI.Ipass
COMMON_BLOCK struct {
long E[96];
long F[96];
long Ix[1001];
long Ir[1001];
long Zone[300];
} ISPACE ;
#define e ISPACE.E
#define f ISPACE.F
#define ix ISPACE.Ix
#define ir ISPACE.Ir
#define zone ISPACE.Zone
COMMON_BLOCK struct {
double U[1001];
double V[1001];
double W[1001];
double X[1001];
double Y[1001];
double Z[1001];
double G[1001];
double Du1[101];
double Du2[101];
double Du3[101];
double Grd[1001];
double Dex[1001];
double Xi[1001];
double Ex[1001];
double Ex1[1001];
double Dex1[1001];
double Vx[1001];
double Xx[1001];
double Rx[1001];
double Rh[2048];
double Vsp[101];
double Vstp[101];
double Vxne[101];
double Vxnd[101];
double Ve3[101];
double Vlr[101];
double Vlin[101];
double B5[101];
double Plan[300];
double D[300];
double Sa[101];
double Sb[101];
} SPACE1 ;
#define u SPACE1.U
#define v SPACE1.V
#define w SPACE1.W
#define x SPACE1.X
#define y SPACE1.Y
#define z SPACE1.Z
#define g SPACE1.G
#define du1 SPACE1.Du1
#define du2 SPACE1.Du2
#define du3 SPACE1.Du3
#define grd SPACE1.Grd
#define dex SPACE1.Dex
#define xi SPACE1.Xi
#define ex SPACE1.Ex
#define ex1 SPACE1.Ex1
#define dex1 SPACE1.Dex1
#define vx SPACE1.Vx
#define xx SPACE1.Xx
#define rx SPACE1.Rx
#define rh SPACE1.Rh
#define vsp SPACE1.Vsp
#define vstp SPACE1.Vstp
#define vxne SPACE1.Vxne
#define vxnd SPACE1.Vxnd
#define ve3 SPACE1.Ve3
#define vlr SPACE1.Vlr
#define vlin SPACE1.Vlin
#define b5 SPACE1.B5
#define plan SPACE1.Plan
#define d SPACE1.D
#define sa SPACE1.Sa
#define sb SPACE1.Sb
COMMON_BLOCK struct {
double P[512][4];
double Px[101][25];
double Cx[101][25];
double Vy[25][101];
double Vh[7][101];
double Vf[7][101];
double Vg[7][101];
double Vs[7][101];
double Za[7][101];
double Zp[7][101];
double Zq[7][101];
double Zr[7][101];
double Zm[7][101];
double Zb[7][101];
double Zu[7][101];
double Zv[7][101];
double Zz[7][101];
double B[64][64];
double C[64][64];
double H[64][64];
double U1[2][101][5];
double U2[2][101][5];
double U3[2][101][5];
} SPACE2 ;
#define p SPACE2.P
#define px SPACE2.Px
#define cx SPACE2.Cx
#define vy SPACE2.Vy
#define vh SPACE2.Vh
#define vf SPACE2.Vf
#define vg SPACE2.Vg
#define vs SPACE2.Vs
#define za SPACE2.Za
#define zp SPACE2.Zp
#define zq SPACE2.Zq
#define zr SPACE2.Zr
#define zm SPACE2.Zm
#define zb SPACE2.Zb
#define zu SPACE2.Zu
#define zv SPACE2.Zv
#define zz SPACE2.Zz
#define b SPACE2.B
#define c SPACE2.C
#define h SPACE2.H
#define u1 SPACE2.U1
#define u2 SPACE2.U2
#define u3 SPACE2.U3
/* KERNEL routine */
CALLED_FROM_FORTRAN KERNEL( TK )
double TK[6];
{
#pragma nodyneqv
#pragma o=i
long argument , k , l , ipnt , ipntp , i;
long lw , j , nl1 , nl2 , kx , ky , ip , kn;
long i1 , j1 , i2 , j2 , nz , ink , jn , kb5i;
long ii , lb , j4 , ng;
double tmp , temp, sum, som;
char name[8];
/*
* Prologue
*/
for(k=0; k<7; k++) name[k] = "kernel"[k];
TRACE( &name );
SPACE();
argument = 0;
TEST( &argument );
/*
*******************************************************************
* Kernel 1 -- hydro fragment
*******************************************************************
* DO 1 L = 1,Loop
* DO 1 k = 1,n
* 1 X(k)= Q + Y(k)*(R*ZX(k+10) + T*ZX(k+11))
*/
for ( l=1 ; l<=loop ; l++ ) {
for ( k=0 ; k<n ; k++ ) {
x[k] = q + y[k]*( r*z[k+10] + t*z[k+11] );
}
}
argument = 1;
TEST( &argument );
/*
*******************************************************************
* Kernel 2 -- ICCG excerpt (Incomplete Cholesky Conjugate Gradient)
*******************************************************************
* DO 200 L= 1,Loop
* II= n
* IPNTP= 0
*222 IPNT= IPNTP
* IPNTP= IPNTP+II
* II= II/2
* i= IPNTP
CDIR$ IVDEP
* DO 2 k= IPNT+2,IPNTP,2
* i= i+1
* 2 X(i)= X(k) - V(k)*X(k-1) - V(k+1)*X(k+1)
* IF( II.GT.1) GO TO 222
*200 CONTINUE
*/
for ( l=1 ; l<=loop ; l++ ) {
ii = n;
ipntp = 0;
do {
ipnt = ipntp;
ipntp += ii;
ii /= 2;
i = ipntp - 1;
#pragma nohazard
for ( k=ipnt+1 ; k<ipntp ; k=k+2 ) {
i++;
x[i] = x[k] - v[k ]*x[k-1] - v[k+1]*x[k+1];
}
} while ( ii>0 );
}
argument = 2;
TEST( &argument );
/*
*******************************************************************
* Kernel 3 -- inner product
*******************************************************************
* DO 3 L= 1,Loop
* Q= 0.0
* DO 3 k= 1,n
* 3 Q= Q + Z(k)*X(k)
*/
for ( l=1 ; l<=loop ; l++ ) {
q = 0.0;
for ( k=0 ; k<n ; k++ ) {
q += z[k]*x[k];
}
}
argument = 3;
TEST( &argument );
/*
*******************************************************************
* Kernel 4 -- banded linear equations
*******************************************************************
* m= (1001-7)/2
* DO 444 L= 1,Loop
* DO 444 k= 7,1001,m
* lw= k-6
* temp= X(k-1)
CDIR$ IVDEP
* DO 4 j= 5,n,5
* temp = temp - XZ(lw)*Y(j)
* 4 lw= lw+1
* X(k-1)= Y(5)*temp
*444 CONTINUE
*/
m = ( 1001-7 )/2;
for ( l=1 ; l<=loop ; l++ ) {
for ( k=6 ; k<1001 ; k=k+m ) {
lw = k - 6;
temp = x[k-1];
#pragma nohazard
for ( j=4 ; j<n ; j=j+5 ) {
temp -= x[lw]*y[j];
lw++;
}
x[k-1] = y[4]*temp;
}
}
argument = 4;
TEST( &argument );
/*
*******************************************************************
* Kernel 5 -- tri-diagonal elimination, below diagonal
*******************************************************************
* DO 5 L = 1,Loop
* DO 5 i = 2,n
* 5 X(i)= Z(i)*(Y(i) - X(i-1))
*/
for ( l=1 ; l<=loop ; l++ ) {
for ( i=1 ; i<n ; i++ ) {
x[i] = z[i]*( y[i] - x[i-1] );
}
}
argument = 5;
TEST( &argument );
/*
*******************************************************************
* Kernel 6 -- general linear recurrence equations
*******************************************************************
* DO 6 L= 1,Loop
* DO 6 i= 2,n
* DO 6 k= 1,i-1
* W(i)= W(i) + B(i,k) * W(i-k)
* 6 CONTINUE
*/
for ( l=1 ; l<=loop ; l++ ) {
for ( i=1 ; i<n ; i++ ) {
for ( k=0 ; k<i ; k++ ) {
w[i] += b[k][i] * w[(i-k)-1];
}
}
}
argument = 6;
TEST( &argument );
/*
*******************************************************************
* Kernel 7 -- equation of state fragment
*******************************************************************
* DO 7 L= 1,Loop
* DO 7 k= 1,n
* X(k)= U(k ) + R*( Z(k ) + R*Y(k )) +
* . T*( U(k+3) + R*( U(k+2) + R*U(k+1)) +
* . T*( U(k+6) + R*( U(k+5) + R*U(k+4))))
* 7 CONTINUE
*/
for ( l=1 ; l<=loop ; l++ ) {
#pragma nohazard
for ( k=0 ; k<n ; k++ ) {
x[k] = u[k] + r*( z[k] + r*y[k] ) +
t*( u[k+3] + r*( u[k+2] + r*u[k+1] ) +
t*( u[k+6] + r*( u[k+5] + r*u[k+4] ) ) );
}
}
argument = 7;
TEST( &argument );
/*
*******************************************************************
* Kernel 8 -- ADI integration
*******************************************************************
* DO 8 L = 1,Loop
* nl1 = 1
* nl2 = 2
* DO 8 kx = 2,3
CDIR$ IVDEP
* DO 8 ky = 2,n
* DU1(ky)=U1(kx,ky+1,nl1) - U1(kx,ky-1,nl1)
* DU2(ky)=U2(kx,ky+1,nl1) - U2(kx,ky-1,nl1)
* DU3(ky)=U3(kx,ky+1,nl1) - U3(kx,ky-1,nl1)
* U1(kx,ky,nl2)=U1(kx,ky,nl1) +A11*DU1(ky) +A12*DU2(ky) +A13*DU3(ky)
* . + SIG*(U1(kx+1,ky,nl1) -2.*U1(kx,ky,nl1) +U1(kx-1,ky,nl1))
* U2(kx,ky,nl2)=U2(kx,ky,nl1) +A21*DU1(ky) +A22*DU2(ky) +A23*DU3(ky)
* . + SIG*(U2(kx+1,ky,nl1) -2.*U2(kx,ky,nl1) +U2(kx-1,ky,nl1))
* U3(kx,ky,nl2)=U3(kx,ky,nl1) +A31*DU1(ky) +A32*DU2(ky) +A33*DU3(ky)
* . + SIG*(U3(kx+1,ky,nl1) -2.*U3(kx,ky,nl1) +U3(kx-1,ky,nl1))
* 8 CONTINUE
*/
for ( l=1 ; l<=loop ; l++ ) {
nl1 = 0;
nl2 = 1;
for ( kx=1 ; kx<3 ; kx++ ){
#pragma nohazard
for ( ky=1 ; ky<n ; ky++ ) {
du1[ky] = u1[nl1][ky+1][kx] - u1[nl1][ky-1][kx];
du2[ky] = u2[nl1][ky+1][kx] - u2[nl1][ky-1][kx];
du3[ky] = u3[nl1][ky+1][kx] - u3[nl1][ky-1][kx];
u1[nl2][ky][kx]=
u1[nl1][ky][kx]+a11*du1[ky]+a12*du2[ky]+a13*du3[ky] + sig*
(u1[nl1][ky][kx+1]-2.0*u1[nl1][ky][kx]+u1[nl1][ky][kx-1]);
u2[nl2][ky][kx]=
u2[nl1][ky][kx]+a21*du1[ky]+a22*du2[ky]+a23*du3[ky] + sig*
(u2[nl1][ky][kx+1]-2.0*u2[nl1][ky][kx]+u2[nl1][ky][kx-1]);
u3[nl2][ky][kx]=
u3[nl1][ky][kx]+a31*du1[ky]+a32*du2[ky]+a33*du3[ky] + sig*
(u3[nl1][ky][kx+1]-2.0*u3[nl1][ky][kx]+u3[nl1][ky][kx-1]);
}
}
}
argument = 8;
TEST( &argument );
/*
*******************************************************************
* Kernel 9 -- integrate predictors
*******************************************************************
* DO 9 L = 1,Loop
* DO 9 i = 1,n
* PX( 1,i)= DM28*PX(13,i) + DM27*PX(12,i) + DM26*PX(11,i) +
* . DM25*PX(10,i) + DM24*PX( 9,i) + DM23*PX( 8,i) +
* . DM22*PX( 7,i) + C0*(PX( 5,i) + PX( 6,i))+ PX( 3,i)
* 9 CONTINUE
*/
for ( l=1 ; l<=loop ; l++ ) {
for ( i=0 ; i<n ; i++ ) {
px[i][0] = dm28*px[i][12] + dm27*px[i][11] + dm26*px[i][10] +
dm25*px[i][ 9] + dm24*px[i][ 8] + dm23*px[i][ 7] +
dm22*px[i][ 6] + c0*( px[i][ 4] + px[i][ 5]) + px[i][ 2];
}
}
argument = 9;
TEST( &argument );
/*
*******************************************************************
* Kernel 10 -- difference predictors
*******************************************************************
* DO 10 L= 1,Loop
* DO 10 i= 1,n
* AR = CX(5,i)
* BR = AR - PX(5,i)
* PX(5,i) = AR
* CR = BR - PX(6,i)
* PX(6,i) = BR
* AR = CR - PX(7,i)
* PX(7,i) = CR
* BR = AR - PX(8,i)
* PX(8,i) = AR
* CR = BR - PX(9,i)
* PX(9,i) = BR
* AR = CR - PX(10,i)
* PX(10,i)= CR
* BR = AR - PX(11,i)
* PX(11,i)= AR
* CR = BR - PX(12,i)
* PX(12,i)= BR
* PX(14,i)= CR - PX(13,i)
* PX(13,i)= CR
* 10 CONTINUE
*/
for ( l=1 ; l<=loop ; l++ ) {
for ( i=0 ; i<n ; i++ ) {
ar = cx[i][ 4];
br = ar - px[i][ 4];
px[i][ 4] = ar;
cr = br - px[i][ 5];
px[i][ 5] = br;
ar = cr - px[i][ 6];
px[i][ 6] = cr;
br = ar - px[i][ 7];
px[i][ 7] = ar;
cr = br - px[i][ 8];
px[i][ 8] = br;
ar = cr - px[i][ 9];
px[i][ 9] = cr;
br = ar - px[i][10];
px[i][10] = ar;
cr = br - px[i][11];
px[i][11] = br;
px[i][13] = cr - px[i][12];
px[i][12] = cr;
}
}
argument = 10;
TEST( &argument );
/*
*******************************************************************
* Kernel 11 -- first sum
*******************************************************************
* DO 11 L = 1,Loop
* X(1)= Y(1)
* DO 11 k = 2,n
* 11 X(k)= X(k-1) + Y(k)
*/
for ( l=1 ; l<=loop ; l++ ) {
x[0] = y[0];
for ( k=1 ; k<n ; k++ ) {
x[k] = x[k-1] + y[k];
}
}
argument = 11;
TEST( &argument );
/*
*******************************************************************
* Kernel 12 -- first difference
*******************************************************************
* DO 12 L = 1,Loop
* DO 12 k = 1,n
* 12 X(k)= Y(k+1) - Y(k)
*/
for ( l=1 ; l<=loop ; l++ ) {
for ( k=0 ; k<n ; k++ ) {
x[k] = y[k+1] - y[k];
}
}
argument = 12;
TEST( &argument );
/*
*******************************************************************
* Kernel 13 -- 2-D PIC (Particle In Cell)
*******************************************************************
* DO 13 L= 1,Loop
* DO 13 ip= 1,n
* i1= P(1,ip)
* j1= P(2,ip)
* i1= 1 + MOD2N(i1,64)
* j1= 1 + MOD2N(j1,64)
* P(3,ip)= P(3,ip) + B(i1,j1)
* P(4,ip)= P(4,ip) + C(i1,j1)
* P(1,ip)= P(1,ip) + P(3,ip)
* P(2,ip)= P(2,ip) + P(4,ip)
* i2= P(1,ip)
* j2= P(2,ip)
* i2= MOD2N(i2,64)
* j2= MOD2N(j2,64)
* P(1,ip)= P(1,ip) + Y(i2+32)
* P(2,ip)= P(2,ip) + Z(j2+32)
* i2= i2 + E(i2+32)
* j2= j2 + F(j2+32)
* H(i2,j2)= H(i2,j2) + 1.0
* 13 CONTINUE
*/
for ( l=1 ; l<=loop ; l++ ) {
for ( ip=0 ; ip<n ; ip++ ) {
i1 = p[ip][0];
j1 = p[ip][1];
i1 &= 64-1;
j1 &= 64-1;
p[ip][2] += b[j1][i1];
p[ip][3] += c[j1][i1];
p[ip][0] += p[ip][2];
p[ip][1] += p[ip][3];
i2 = p[ip][0];
j2 = p[ip][1];
i2 = ( i2 & 64-1 ) - 1 ;
j2 = ( j2 & 64-1 ) - 1 ;
p[ip][0] += y[i2+32];
p[ip][1] += z[j2+32];
i2 += e[i2+32];
j2 += f[j2+32];
h[j2][i2] += 1.0;
}
}
argument = 13;
TEST( &argument );
/*
*******************************************************************
* Kernel 14 -- 1-D PIC (Particle In Cell)
*******************************************************************
* DO 14 L= 1,Loop
* DO 141 k= 1,n
* VX(k)= 0.0
* XX(k)= 0.0
* IX(k)= INT( GRD(k))
* XI(k)= REAL( IX(k))
* EX1(k)= EX ( IX(k))
* DEX1(k)= DEX ( IX(k))
*41 CONTINUE
* DO 142 k= 1,n
* VX(k)= VX(k) + EX1(k) + (XX(k) - XI(k))*DEX1(k)
* XX(k)= XX(k) + VX(k) + FLX
* IR(k)= XX(k)
* RX(k)= XX(k) - IR(k)
* IR(k)= MOD2N( IR(k),2048) + 1
* XX(k)= RX(k) + IR(k)
*42 CONTINUE
* DO 14 k= 1,n
* RH(IR(k) )= RH(IR(k) ) + 1.0 - RX(k)
* RH(IR(k)+1)= RH(IR(k)+1) + RX(k)
*14 CONTINUE
*/
for ( l=1 ; l<=loop ; l++ ) {
for ( k=0 ; k<n ; k++ ) {
vx[k] = 0.0;
xx[k] = 0.0;
ix[k] = (long) grd[k];
xi[k] = (double) ix[k];
ex1[k] = ex[ ix[k] - 1 ];
dex1[k] = dex[ ix[k] - 1 ];
}
for ( k=0 ; k<n ; k++ ) {
vx[k] = vx[k] + ex1[k] + ( xx[k] - xi[k] )*dex1[k];
xx[k] = xx[k] + vx[k] + flx;
ir[k] = xx[k];
rx[k] = xx[k] - ir[k];
ir[k] = ( ir[k] & 2048-1 ) + 1;
xx[k] = rx[k] + ir[k];
}
for ( k=0 ; k<n ; k++ ) {
rh[ ir[k]-1 ] += 1.0 - rx[k];
rh[ ir[k] ] += rx[k];
}
}
argument = 14;
TEST( &argument );
/*
*******************************************************************
* Kernel 15 -- Casual Fortran. Development version
*******************************************************************
* DO 45 L = 1,Loop
* NG= 7
* NZ= n
* AR= 0.053
* BR= 0.073
* 15 DO 45 j = 2,NG
* DO 45 k = 2,NZ
* IF( j-NG) 31,30,30
* 30 VY(k,j)= 0.0
* GO TO 45
* 31 IF( VH(k,j+1) -VH(k,j)) 33,33,32
* 32 T= AR
* GO TO 34
* 33 T= BR
* 34 IF( VF(k,j) -VF(k-1,j)) 35,36,36
* 35 R= MAX( VH(k-1,j), VH(k-1,j+1))
* S= VF(k-1,j)
* GO TO 37
* 36 R= MAX( VH(k,j), VH(k,j+1))
* S= VF(k,j)
* 37 VY(k,j)= SQRT( VG(k,j)**2 +R*R)*T/S
* 38 IF( k-NZ) 40,39,39
* 39 VS(k,j)= 0.
* GO TO 45
* 40 IF( VF(k,j) -VF(k,j-1)) 41,42,42
* 41 R= MAX( VG(k,j-1), VG(k+1,j-1))
* S= VF(k,j-1)
* T= BR
* GO TO 43
* 42 R= MAX( VG(k,j), VG(k+1,j))
* S= VF(k,j)
* T= AR
* 43 VS(k,j)= SQRT( VH(k,j)**2 +R*R)*T/S
* 45 CONTINUE
*/
#pragma intrinsic sqrt
for ( l=1 ; l<=loop ; l++ ) {
ng = 7;
nz = n;
ar = 0.053;
br = 0.073;
for ( j=1 ; j<ng ; j++ ) {
for ( k=1 ; k<nz ; k++ ) {
if ( (j+1) >= ng ) {
vy[j][k] = 0.0;
continue;
}
if ( vh[j+1][k] > vh[j][k] ) {
t = ar;
}
else {
t = br;
}
if ( vf[j][k] < vf[j][k-1] ) {
if ( vh[j][k-1] > vh[j+1][k-1] )
r = vh[j][k-1];
else
r = vh[j+1][k-1];
s = vf[j][k-1];
}
else {
if ( vh[j][k] > vh[j+1][k] )
r = vh[j][k];
else
r = vh[j+1][k];
s = vf[j][k];
}
vy[j][k] = sqrt( vg[j][k]*vg[j][k] + r*r )* t/s;
if ( (k+1) >= nz ) {
vs[j][k] = 0.0;
continue;
}
if ( vf[j][k] < vf[j-1][k] ) {
if ( vg[j-1][k] > vg[j-1][k+1] )
r = vg[j-1][k];
else
r = vg[j-1][k+1];
s = vf[j-1][k];
t = br;
}
else {
if ( vg[j][k] > vg[j][k+1] )
r = vg[j][k];
else
r = vg[j][k+1];
s = vf[j][k];
t = ar;
}
vs[j][k] = sqrt( vh[j][k]*vh[j][k] + r*r )* t / s;
}
}
}
argument = 15;
TEST( &argument );
/*
*******************************************************************
* Kernel 16 -- Monte Carlo search loop
*******************************************************************
* II= n/3
* LB= II+II
* k2= 0
* k3= 0
* DO 485 L= 1,Loop
* m= 1
*405 i1= m
*410 j2= (n+n)*(m-1)+1
* DO 470 k= 1,n
* k2= k2+1
* j4= j2+k+k
* j5= ZONE(j4)
* IF( j5-n ) 420,475,450
*415 IF( j5-n+II ) 430,425,425
*420 IF( j5-n+LB ) 435,415,415
*425 IF( PLAN(j5)-R) 445,480,440
*430 IF( PLAN(j5)-S) 445,480,440
*435 IF( PLAN(j5)-T) 445,480,440
*440 IF( ZONE(j4-1)) 455,485,470
*445 IF( ZONE(j4-1)) 470,485,455
*450 k3= k3+1
* IF( D(j5)-(D(j5-1)*(T-D(j5-2))**2+(S-D(j5-3))**2
* . +(R-D(j5-4))**2)) 445,480,440
*455 m= m+1
* IF( m-ZONE(1) ) 465,465,460
*460 m= 1
*465 IF( i1-m) 410,480,410
*470 CONTINUE
*475 CONTINUE
*480 CONTINUE
*485 CONTINUE
*/
ii = n / 3;
lb = ii + ii;
k3 = k2 = 0;
for ( l=1 ; l<=loop ; l++ ) {
i1 = m = 1;
label410:
j2 = ( n + n )*( m - 1 ) + 1;
for ( k=1 ; k<=n ; k++ ) {
k2++;
j4 = j2 + k + k;
j5 = zone[j4-1];
if ( j5 < n ) {
if ( j5+lb < n ) { /* 420 */
tmp = plan[j5-1] - t; /* 435 */
} else {
if ( j5+ii < n ) { /* 415 */
tmp = plan[j5-1] - s; /* 430 */
} else {
tmp = plan[j5-1] - r; /* 425 */
}
}
} else if( j5 == n ) {
break; /* 475 */
} else {
k3++; /* 450 */
tmp=(d[j5-1]-(d[j5-2]*(t-d[j5-3])*(t-d[j5-3])+(s-d[j5-4])*
(s-d[j5-4])+(r-d[j5-5])*(r-d[j5-5])));
}
if ( tmp < 0.0 ) {
if ( zone[j4-2] < 0 ) /* 445 */
continue; /* 470 */
else if ( !zone[j4-2] )
break; /* 480 */
} else if ( tmp ) {
if ( zone[j4-2] > 0 ) /* 440 */
continue; /* 470 */
else if ( !zone[j4-2] )
break; /* 480 */
} else break; /* 485 */
m++; /* 455 */
if ( m > zone[0] )
m = 1; /* 460 */
if ( i1-m ) /* 465 */
goto label410;
else
break;
}
}
argument = 16;
TEST( &argument );
/*
*******************************************************************
* Kernel 17 -- implicit, conditional computation
*******************************************************************
* DO 62 L= 1,Loop
* i= n
* j= 1
* INK= -1
* SCALE= 5./3.
* XNM= 1./3.
* E6= 1.03/3.07
* GO TO 61
*60 E6= XNM*VSP(i)+VSTP(i)
* VXNE(i)= E6
* XNM= E6
* VE3(i)= E6
* i= i+INK
* IF( i.EQ.j) GO TO 62
*61 E3= XNM*VLR(i) +VLIN(i)
* XNEI= VXNE(i)
* VXND(i)= E6
* XNC= SCALE*E3
* IF( XNM .GT.XNC) GO TO 60
* IF( XNEI.GT.XNC) GO TO 60
* VE3(i)= E3
* E6= E3+E3-XNM
* VXNE(i)= E3+E3-XNEI
* XNM= E6
* i= i+INK
* IF( i.NE.j) GO TO 61
* 62 CONTINUE
*/
for ( l=1 ; l<=loop ; l++ ) {
i = n-1;
j = 0;
ink = -1;
scale = 5.0 / 3.0;
xnm = 1.0 / 3.0;
e6 = 1.03 / 3.07;
goto l61;
l60: e6 = xnm*vsp[i] + vstp[i];
vxne[i] = e6;
xnm = e6;
ve3[i] = e6;
i += ink;
if ( i==j ) goto l62;
l61: e3 = xnm*vlr[i] + vlin[i];
xnei = vxne[i];
vxnd[i] = e6;
xnc = scale*e3;
if ( xnm > xnc ) goto l60;
if ( xnei > xnc ) goto l60;
ve3[i] = e3;
e6 = e3 + e3 - xnm;
vxne[i] = e3 + e3 - xnei;
xnm = e6;
i += ink;
if ( i != j ) goto l61;
l62:;
}
argument = 17;
TEST( &argument );
/*
*******************************************************************
* Kernel 18 - 2-D explicit hydrodynamics fragment
*******************************************************************
* DO 75 L= 1,Loop
* T= 0.0037
* S= 0.0041
* KN= 6
* JN= n
* DO 70 k= 2,KN
* DO 70 j= 2,JN
* ZA(j,k)= (ZP(j-1,k+1)+ZQ(j-1,k+1)-ZP(j-1,k)-ZQ(j-1,k))
* . *(ZR(j,k)+ZR(j-1,k))/(ZM(j-1,k)+ZM(j-1,k+1))
* ZB(j,k)= (ZP(j-1,k)+ZQ(j-1,k)-ZP(j,k)-ZQ(j,k))
* . *(ZR(j,k)+ZR(j,k-1))/(ZM(j,k)+ZM(j-1,k))
* 70 CONTINUE
* DO 72 k= 2,KN
* DO 72 j= 2,JN
* ZU(j,k)= ZU(j,k)+S*(ZA(j,k)*(ZZ(j,k)-ZZ(j+1,k))
* . -ZA(j-1,k) *(ZZ(j,k)-ZZ(j-1,k))
* . -ZB(j,k) *(ZZ(j,k)-ZZ(j,k-1))
* . +ZB(j,k+1) *(ZZ(j,k)-ZZ(j,k+1)))
* ZV(j,k)= ZV(j,k)+S*(ZA(j,k)*(ZR(j,k)-ZR(j+1,k))
* . -ZA(j-1,k) *(ZR(j,k)-ZR(j-1,k))
* . -ZB(j,k) *(ZR(j,k)-ZR(j,k-1))
* . +ZB(j,k+1) *(ZR(j,k)-ZR(j,k+1)))
* 72 CONTINUE
* DO 75 k= 2,KN
* DO 75 j= 2,JN
* ZR(j,k)= ZR(j,k)+T*ZU(j,k)
* ZZ(j,k)= ZZ(j,k)+T*ZV(j,k)
* 75 CONTINUE
*/
for ( l=1 ; l<=loop ; l++ ) {
t = 0.0037;
s = 0.0041;
kn = 6;
jn = n;
for ( k=1 ; k<kn ; k++ ) {
#pragma nohazard
for ( j=1 ; j<jn ; j++ ) {
za[k][j] = ( zp[k+1][j-1] +zq[k+1][j-1] -zp[k][j-1] -zq[k][j-1] )*
( zr[k][j] +zr[k][j-1] ) / ( zm[k][j-1] +zm[k+1][j-1]);
zb[k][j] = ( zp[k][j-1] +zq[k][j-1] -zp[k][j] -zq[k][j] ) *
( zr[k][j] +zr[k-1][j] ) / ( zm[k][j] +zm[k][j-1]);
}
}
for ( k=1 ; k<kn ; k++ ) {
#pragma nohazard
for ( j=1 ; j<jn ; j++ ) {
zu[k][j] += s*( za[k][j] *( zz[k][j] - zz[k][j+1] ) -
za[k][j-1] *( zz[k][j] - zz[k][j-1] ) -
zb[k][j] *( zz[k][j] - zz[k-1][j] ) +
zb[k+1][j] *( zz[k][j] - zz[k+1][j] ) );
zv[k][j] += s*( za[k][j] *( zr[k][j] - zr[k][j+1] ) -
za[k][j-1] *( zr[k][j] - zr[k][j-1] ) -
zb[k][j] *( zr[k][j] - zr[k-1][j] ) +
zb[k+1][j] *( zr[k][j] - zr[k+1][j] ) );
}
}
for ( k=1 ; k<kn ; k++ ) {
#pragma nohazard
for ( j=1 ; j<jn ; j++ ) {
zr[k][j] = zr[k][j] + t*zu[k][j];
zz[k][j] = zz[k][j] + t*zv[k][j];
}
}
}
argument = 18;
TEST( &argument );
/*
*******************************************************************
* Kernel 19 -- general linear recurrence equations
*******************************************************************
* KB5I= 0
* DO 194 L= 1,Loop
* DO 191 k= 1,n
* B5(k+KB5I)= SA(k) +STB5*SB(k)
* STB5= B5(k+KB5I) -STB5
*191 CONTINUE
*192 DO 193 i= 1,n
* k= n-i+1
* B5(k+KB5I)= SA(k) +STB5*SB(k)
* STB5= B5(k+KB5I) -STB5
*193 CONTINUE
*194 CONTINUE
*/
kb5i = 0;
for ( l=1 ; l<=loop ; l++ ) {
for ( k=0 ; k<n ; k++ ) {
b5[k+kb5i] = sa[k] + stb5*sb[k];
stb5 = b5[k+kb5i] - stb5;
}
for ( i=1 ; i<=n ; i++ ) {
k = n - i ;
b5[k+kb5i] = sa[k] + stb5*sb[k];
stb5 = b5[k+kb5i] - stb5;
}
}
argument = 19;
TEST( &argument );
/*
*******************************************************************
* Kernel 20 -- Discrete ordinates transport, conditional recurrence on xx
*******************************************************************
* DO 20 L= 1,Loop
* DO 20 k= 1,n
* DI= Y(k)-G(k)/( XX(k)+DK)
* DN= 0.2
* IF( DI.NE.0.0) DN= MAX( S,MIN( Z(k)/DI, T))
* X(k)= ((W(k)+V(k)*DN)* XX(k)+U(k))/(VX(k)+V(k)*DN)
* XX(k+1)= (X(k)- XX(k))*DN+ XX(k)
* 20 CONTINUE
*/
for ( l=1 ; l<=loop ; l++ ) {
for ( k=0 ; k<n ; k++ ) {
di = y[k] - g[k] / ( xx[k] + dk );
dn = 0.2;
if ( di ) {
dn = z[k]/di ;
if ( t < dn ) dn = t;
if ( s > dn ) dn = s;
}
x[k] = ( ( w[k] + v[k]*dn )* xx[k] + u[k] ) / ( vx[k] + v[k]*dn );
xx[k+1] = ( x[k] - xx[k] )* dn + xx[k];
}
}
argument = 20;
TEST( &argument );
/*
*******************************************************************
* Kernel 21 -- matrix*matrix product
*******************************************************************
* DO 21 L= 1,Loop
* DO 21 k= 1,25
* DO 21 i= 1,25
* DO 21 j= 1,n
* PX(i,j)= PX(i,j) +VY(i,k) * CX(k,j)
* 21 CONTINUE
*/
for ( l=1 ; l<=loop ; l++ ) {
for ( k=0 ; k<25 ; k++ ) {
for ( i=0 ; i<25 ; i++ ) {
#pragma nohazard
for ( j=0 ; j<n ; j++ ) {
px[j][i] += vy[k][i] * cx[j][k];
}
}
}
}
argument = 21;
TEST( &argument );
/*
*******************************************************************
* Kernel 22 -- Planckian distribution
*******************************************************************
* EXPMAX= 20.0
* U(n)= 0.99*EXPMAX*V(n)
* DO 22 L= 1,Loop
* DO 22 k= 1,n
* Y(k)= U(k)/V(k)
* W(k)= X(k)/( EXP( Y(k)) -1.0)
* 22 CONTINUE
*/
#pragma intrinsic exp
expmax = 20.0;
u[n-1] = 0.99*expmax*v[n-1];
for ( l=1 ; l<=loop ; l++ ) {
for ( k=0 ; k<n ; k++ ) {
y[k] = u[k] / v[k];
w[k] = x[k] / ( exp( y[k] ) -1.0 );
}
}
argument = 22;
TEST( &argument );
/*
*******************************************************************
* Kernel 23 -- 2-D implicit hydrodynamics fragment
*******************************************************************
* DO 23 L= 1,Loop
* DO 23 j= 2,6
* DO 23 k= 2,n
* QA= ZA(k,j+1)*ZR(k,j) +ZA(k,j-1)*ZB(k,j) +
* . ZA(k+1,j)*ZU(k,j) +ZA(k-1,j)*ZV(k,j) +ZZ(k,j)
* 23 ZA(k,j)= ZA(k,j) +.175*(QA -ZA(k,j))
*/
for ( l=1 ; l<=loop ; l++ ) {
for ( j=1 ; j<6 ; j++ ) {
for ( k=1 ; k<n ; k++ ) {
qa = za[j+1][k]*zr[j][k] + za[j-1][k]*zb[j][k] +
za[j][k+1]*zu[j][k] + za[j][k-1]*zv[j][k] + zz[j][k];
za[j][k] += 0.175*( qa - za[j][k] );
}
}
}
argument = 23;
TEST( &argument );
/*
*******************************************************************
* Kernel 24 -- find location of first minimum in array
*******************************************************************
* X( n/2)= -1.0E+10
* DO 24 L= 1,Loop
* m= 1
* DO 24 k= 2,n
* IF( X(k).LT.X(m)) m= k
* 24 CONTINUE
*/
x[n/2] = -1.0e+10;
for ( l=1 ; l<=loop ; l++ ) {
m = 0;
for ( k=1 ; k<n ; k++ ) {
if ( x[k] < x[m] ) m = k;
}
}
argument = 24;
TEST( &argument );
/*
* Epilogue
*/
for ( k=0 ; k<mk ; k++ ) {
times[k][il-1][jr-1] = time[k];
fopn[k][il-1][jr-1] = flopn[k];
terrs[k][il-1][jr-1] = terr1[k];
npfs[k][il-1][jr-1] = npfs1[k];
csums[k][il-1][jr-1] = csum[k];
dos[k][il-1][jr-1] = total[k];
}
sum = 0.0;
for ( k=0 ; k<mk ; k++ ) {
sum += time[k];
}
TK[0] += sum;
som = 0.0;
for ( k=0 ; k<mk ; k++ ) {
som += flopn[k]*total[k];
}
TK[1] += som;
TRACK( &name );
return;
}
/* End of File */
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment