📄 livermore.c
字号:
/*
* 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
*/
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -