WOW! - amazing stuff!...SO impressive! - do tell us some more of your exploits! ( oh!..looking back - I see you already have! ) 🙄
Always happy to accommodate.
PROGRAM FSP3D
C HOST CONTROL PROGRAM FOR CLOUD MODEL.
C THREE-DIMENSIONAL FSP PROCESS:
C 3-D PULSE IS A STRAIGHT CYLINDER WHOSE BASE IS A SPHERE.
C IN THIS CASE, THE DISTRIBUTION OF THE PULSE AREAS SHOULD
C BE PR(A>A) = 1/A FOR A>1. THE CENTERS OF THE CIRCLES ARE
C PLACED UNIFORMLY AT RANDOM IN A SQUARE L * L * L.
C MANDELBROT'S FRACTAL SUM OF PULSES ALGORITHM
C WITH EXPONENTIAL DECAY TO GET RID OF PULSE EDGES.
C
C S= DECAY FACTOR: EXP( -( U/ RO)** 2S)
C TOTAL INTERVAL LENGTH = L
C NUMINC = (SMALL L) IN PAPER
C NU PULSE CENTERS/UNIT LENGTH
C NUMCYL= L*NU
C IMPLICIT NONE
INTEGER* 4 NUMCYL, LX, LY, LZ, S, IMGFIL
INTEGER* 4 NUMXINC, NUMYINC, J, L, N, IWATCH, BUFRSZ
INTEGER* 4 IMAGES
REAL* 4 X, R, ALPHA, NU, SUM, LINCX, LINCY, DLX, DLY, THRESH, PI
REAL* 4 MXCLRS, ZINIT, ZDELTA
PARAMETER ( NUMXINC= 480, NUMYINC= 480, BUFRSZ= NUMXINC* 2)
PARAMETER ( NU= 0.025, MXCLRS= 255.0)
PARAMETER ( ALPHA= 5./ 3., NUMCYL= 2000, IMAGES= 6)
PARAMETER ( LX= 4800, LY= 4800, LZ= 4800, PI= 3.1415927)
PARAMETER ( ZINIT= 1200.0, ZDELTA= 400.0)
PARAMETER ( KINIT= 1.025, KDELTA= 0.025)
INTEGER* 4 IMAGNM
REAL* 4 K, SEED, IMGMAX, IMGMIN, ZPOSTN
REAL* 4 CR( NUMXINC, NUMYINC), MSGBFR( BUFRSZ), CRMIN, CRMAX
LOGICAL* 1 IMGBYT( 2, NUMXINC)
INTEGER* 2 IMAGE( NUMXINC)
EQUIVALENCE ( IMAGE( 1), IMGBYT( 1, 1))
CHARACTER* 40 FLNAME
C CUBE SPECIFIC PARAMETERS
INTEGER* 4 INT4, REAL4, REAL8, HOST, ALLNDS
INTEGER* 4 NDMESG, SDMESG, KMESG, ZMESG
PARAMETER ( INT4= 4, REAL4= 4, REAL8= 8)
PARAMETER ( HOST= -32768, ALLNDS= -1)
PARAMETER ( NDMESG= 100, SDMESG= NDMESG, KMESG= NDMESG+ 1)
PARAMETER ( ZMESG= NDMESG+ 2)
INTEGER* 4 PID, CHANNL, PRCSRS
INTEGER* 4 MSGTYP, RECEVD, MSGORG, ORGPID, LENGTH, HSTPID
INTEGER* 4 DTAROW, FINISH, ROMESG, RIMESG, AMESG
PARAMETER ( DTAROW= 1, FINISH= 2, PRCSRS= 4)
PARAMETER ( ROMESG= 3, RIMESG= 4, AMESG= 5)
INTEGER* 4 CUBEDN, ROW, COLUMN, MAXBFR
C EXTERNALS.
INTEGER* 4 MOD, IAND, COPEN
C GET COMM CHANNEL ID, NODE NUMBER, NUMBER OF PROCESSORS.
PID= 0
CHANNL= COPEN( PID)
C GET 'RANDOM' SEED FOR RANDOM NO. GENERATOR.
WRITE( *, 100)
100 FORMAT( ' INPUT A SEED: ')
READ( *, 200) SEED
200 FORMAT( F10.0)
C SEND THE SEED TO ALL NODES.
MSGTYP= SDMESG
CALL SENDMSG( CHANNL, MSGTYP, SEED, REAL4, ALLNDS, PID)
C READ K, SETS THE RING SIZE.
C SEND K TO ALL NODES.
WRITE( *, 101)
101 FORMAT( ' INPUT K (FROM 1 TO 1.4): ')
READ( *, 200) K
MSGTYP= KMESG
CALL SENDMSG( CHANNL, MSGTYP, K, REAL4, ALLNDS, PID)
ZPOSTN= ZINIT
DO 1000 IMAGNM= 1, IMAGES
MSGTYP= ZMESG
CALL SENDMSG( CHANNL, MSGTYP, ZPOSTN, REAL4, ALLNDS, PID)
C PRIME THE MAX AND MIN.
CRMAX= -1.0E30
CRMIN= 1.0E30
CUBEDN= 0
C LISTEN FOR MESSAGES.
MAXBFR= REAL4* BUFRSZ
5 CONTINUE
CALL RECVMSG ( CHANNL, MSGTYP, MSGBFR, MAXBFR, RECEVD,
1 MSGORG, ORGPID)
WRITE( *, *) ' RECEIVED MESSAGE NUMBER ', MSGTYP
IF( MSGTYP .EQ. DTAROW) THEN
C A ROW OF IMAGE DATA HAS BEEN RECEIVED.
C FIND ROW NUMBER.
ROW= INT( MSGBFR( NUMXINC+ 1)+ 0.001)
WRITE( *, *) ' RECEIVED ROW ', ROW, '.'
C TRANSFER TO IMAGE ARRAY.
DO 10 COLUMN= 1, NUMXINC
CR( COLUMN, ROW)= MSGBFR( COLUMN)
10 CONTINUE
ELSE IF( MSGTYP .EQ. FINISH) THEN
C A NODE HAS FINISHED ITS PROCESSING.
CUBEDN= CUBEDN+ 1
C GET MAX AND MIN OF CUBE.
IF( CRMIN .GT. MSGBFR( 1)) CRMIN= MSGBFR( 1)
IF( CRMAX .LT. MSGBFR( 2)) CRMAX= MSGBFR( 2)
ELSE IF( MSGTYP .EQ. ROMESG) THEN
WRITE( *, *) ' RADOUT FOLLOWS:'
WRITE( *, 300)( MSGBFR( COLUMN), COLUMN= 1, 200)
300 FORMAT( 5( 1PG15.5))
ELSE IF( MSGTYP .EQ. RIMESG) THEN
WRITE( *, *) ' RADIN FOLLOWS:'
WRITE( *, 300)( MSGBFR( COLUMN), COLUMN= 1, 200)
ELSE IF( MSGTYP .EQ. AMESG) THEN
WRITE( *, *) ' A FOLLOWS:'
WRITE( *, 300)( MSGBFR( COLUMN), COLUMN= 1, 200)
END IF
C GO BACK FOR MORE MESSAGES IF NODES ARE NOT FINISHED.
IF( CUBEDN .LT. PRCSRS) GO TO 5
C SCALE DATA AND WRITE TO DISK.
WRITE( FLNAME, 400) IMAGNM
400 FORMAT( 'cloud', I2.2, '.img')
OPEN( UNIT= IMGFIL, FILE= FLNAME, FORM= 'UNFORMATTED',
1 ACCESS= 'DIRECT', RECL= NUMXINC)
DO 20 ROW= 1, NUMYINC
DO 15 COLUMN= 1, NUMXINC
IMAGE( COLUMN)= INT( MXCLRS*
1 ( CR( COLUMN, ROW)- CRMIN)/ ( CRMAX- CRMIN))
15 CONTINUE
WRITE( IMGFIL, REC= ROW)
1 ( IMGBYT( 1, COLUMN), COLUMN= 1, NUMXINC)
20 CONTINUE
C WRITE IMAGE SIZE AND IMAGE VARIABLES.
WRITE( IMGFIL, REC= NUMYINC+ 1)
1 NUMXINC, NUMYINC, K, SEED, THRESH, LX, LY, ALPHA, CRMAX, CRMIN,
2 NUMCYL, ZPOSTN
CLOSE( UNIT= IMGFIL)
ZPOSTN= ZPOSTN+ ZDELTA
1000 CONTINUE
C TURN OFF NODES.
K= -1.0
MSGTYP= KMESG
CALL SENDMSG( CHANNL, MSGTYP, K, REAL4, ALLNDS, PID)
STOP 'NORMAL PROGRAM END'
END
C==================================
C NODE CODE.
C LOAD INTO EACH OF FOUR VECTOR PROCESSORS.
PROGRAM FSP2D
C TWO-DIMENSIONAL FSP PROCESS:
C 2-D PULSE IS A STRAIGHT CYLINDER WHOSE BASE IS A CIRCLE.
C IN THIS CASE, THE DISTRIBUTION OF THE PULSE AREAS SHOULD
C BE PR(A>A) = 1/A FOR A>1. THE CENTERS OF THE CIRCLES ARE
C PLACED UNIFORMLY AT RANDOM IN A SQUARE L * L.
C MANDELBROT'S FRACTAL SUM OF PULSES ALGORITHM
C WITH EXPONENTIAL DECAY TO GET RID OF PULSE EDGES.
C
C S= DECAY FACTOR: EXP(-(U/RO)**2S), ALWAYS 1.
C TOTAL INTERVAL LENGTH = L
C NUMINC = (SMALL L) IN PAPER
C NU PULSE CENTERS/UNIT LENGTH
C NUMCYL= L*NU
C IMPLICIT NONE
INTEGER* 4 NUMCYL, LX, LY, S, IMGFIL
INTEGER* 4 NUMXINC, NUMYINC, J, L, N, IWATCH
REAL* 4 LXREAL, LYREAL
REAL* 4 X, R, ALPHA, NU, SUM, LINCX, LINCY, DLX, DLY, THRESH, PI
PARAMETER( NUMXINC= 600, NUMYINC= 400, NU= 0.025)
PARAMETER( ALPHA= 5./ 3., NUMCYL= 20000)
PARAMETER( LX= 4800, LY= 3000, PI= 3.1415927)
PARAMETER ( LXREAL= 4800.0, LYREAL= 3000.0)
C+++ COMMON/ VECSPC/
INTEGER* 4 CRBFSZ
PARAMETER ( CRBFSZ= NUMXINC+ 1)
REAL* 4 DR, INPI, SIGIN, DEL, SIGNBT
REAL* 4 DIS, POSX, POSY, A, K, RADIN, RADOUT
REAL* 4 SEED, IMGMAX, IMGMIN, CR, CRMIN, CRMAX
COMMON/ VECSPC/
1 DR( NUMCYL), SIGIN( NUMCYL),
2 DEL( NUMCYL), POSX( NUMCYL), POSY( NUMCYL), A( NUMCYL),
3 RADIN( NUMCYL), RADOUT( NUMCYL), CR( CRBFSZ),
4 VTEMP1( NUMCYL), VTEMP2( NUMCYL), VTEMP3( NUMCYL),
5 SEED, IMGMAX, IMGMIN, SIGNBT, K, INPI, CRMIN, CRMAX
LOGICAL* 4 LMASK1( NUMCYL), LMASK2( NUMCYL)
INTEGER* 4 ITEMP1( NUMCYL), ITEMP2( NUMCYL)
EQUIVALENCE ( LMASK1( 1), ITEMP1( 1), VTEMP1( 1))
EQUIVALENCE ( LMASK2( 1), ITEMP2( 1), VTEMP2( 1))
C--- END / VECSPC/
C OPTIMIZING VARIABLES.
REAL* 4 STEMP1
C CUBE SPECIFIC VARIABLES.
INTEGER* 4 INT4, REAL4, REAL8, HOST, ALLNDS
INTEGER* 4 NDMESG, SDMESG, KMESG
PARAMETER ( INT4= 4, REAL4= 4, REAL8= 8)
PARAMETER ( HOST= -32768, ALLNDS= -1)
PARAMETER ( NDMESG= 100, SDMESG= NDMESG, KMESG= NDMESG+ 1)
INTEGER* 4 PID, CHANNL, PRCSRS
INTEGER* 4 MSGTYP, RECEVD, MSGORG, ORGPID, LENGTH, HSTPID
INTEGER* 4 XSTRID, YSTRID, ZSTRID
C EXTERNALS.
INTEGER* 4 MOD, IAND
INTEGER* 4 COPEN, CUBEDIM, MYNODE
REAL* 4 SDOT
C PARAMETER (NUMINC=100, NU= .002, L= 200.)
C PARAMETER (ALPHA= 5./3., NUMCYL= 1)
C NODE INITIALIZATION.
C GET COMM CHANNEL ID, NODE NUMBER, NUMBER OF PROCESSORS.
C
PID= 0
CHANNL= COPEN( PID)
NODEID= MYNODE( )
PRCSRS= 2** CUBEDIM( )
C ATTACH SOME SPEED.
CALL ATTACH
C DR IS THE PULSE INTENSITY FOR CYLINDER RADIUS RO.
C FOR A GIVEN RADIUS RO, THE PULSE INTENSITY
C DR IS OF FIXED ABSOLUTE VALUE AND OF RANDOM SIGN.
C GENERATES A RV WITH 1-1/X DISTRIBUTION
C
C GET 'RANDOM' SEED FOR RANDOM NO. GENERATOR FROM THE HOST PROGRAM.
MSGTYP= SDMESG
CALL VPWLEDS( 0)
CALL GREENLED( 1)
CALL RECVW ( CHANNL, MSGTYP, SEED, REAL4, RECEVD, MSGORG, ORGPID)
IF( MSGORG .EQ. HOST) THEN
HSTPID= ORGPID
ELSE
STOP 'BAD MESSAGE'
END IF
C GENERATE ANNULI BY ENTERING THE INNER RADIUS. THE INNER RADIUS
C IS SPECIFIED BY THE VALUE OF K. IF K=1, THEN THE ANNULI ARE
C CIRCLES. K=1.4 RESULTS IN VERY THIN RINGS. RADOUT IS THE RADIUS
C OF THE OUTER CIRCLE, RADIN IS THE RADIUS OF THE INNER CIRCLE.
MSGTYP= KMESG
CALL RECVW ( CHANNL, MSGTYP, K, REAL4, RECEVD, MSGORG, ORGPID)
CALL REDLED( 1)
CALL GREENLED( 0)
C PLACE CYLINDER CENTERS UNIFORMLY AT RANDOM OVER LX X LY RECTANGLE,
C NU PER UNIT AREA, WITH AREA DISTRIBUTED AS 1/RO FOR RO>THRESH
C THE AMPLITUDE DR IS OF FIXED ABSOLUTE VALUE AND OF RANDOM SIGN.
INPI= 1.0/ PI
C GENERATE ALL REQUIRED RANDOM NUMBERS.
CALL VPWLEDS( 1)
CALL FVRAND( SEED, A, NUMCYL)
CALL FVRAND( SEED, POSX, NUMCYL)
CALL FVRAND( SEED, POSY, NUMCYL)
CALL FVRAND( SEED, DR, NUMCYL)
XSTRID= 1
YSTRID= 1
ZSTRID= 1
CALL VPWLEDS( 2)
C BLAST OFF!!!
C DO 10 J= 1, NUMCYL
C A( J)= 1.E5/( 1.0- A( J))
CALL SSSUB( NUMCYL, 1.0, A, XSTRID, A, YSTRID)
CALL SSDIV( NUMCYL, 1.0E5, A, XSTRID, A, YSTRID)
C RADOUT( J)= SQRT( A( J)* INPI)
CALL SSMUL( NUMCYL, INPI, A, XSTRID, RADOUT, YSTRID)
CALL SVSQRT( NUMCYL, RADOUT, XSTRID, RADOUT, YSTRID)
C RADIN( J)= SQRT( A( J)* INPI*( K** 2- 1.0))
STEMP1= INPI* ( K** 2- 1.0)
CALL SSMUL( NUMCYL, STEMP1, A, XSTRID, RADIN, YSTRID)
CALL SVSQRT( NUMCYL, RADIN, XSTRID, RADIN, YSTRID)
C DEL( J)= ( 0.5*( RADIN( J)+ RADOUT( J)))** 2
CALL SSVVPT( NUMCYL, 0.5, RADIN, XSTRID, RADOUT, YSTRID,
1 DEL, ZSTRID)
CALL SVMUL( NUMCYL, DEL, XSTRID, DEL, XSTRID, DEL, XSTRID)
C SIGIN( J)= ( 1.0/( 0.5*( RADOUT( J)- RADIN( J))))** 2
CALL SSVVMT( NUMCYL, 0.5, RADOUT, XSTRID, RADIN, YSTRID,
1 SIGIN, ZSTRID)
CALL SVRECP( NUMCYL, SIGIN, XSTRID, SIGIN, XSTRID)
CALL SVMUL( NUMCYL, SIGIN, XSTRID, SIGIN, XSTRID,
1 SIGIN, XSTRID)
C POSX( J)= POSX( J)* LXREAL
CALL SSMUL( NUMCYL, LXREAL, POSX, XSTRID, POSX, XSTRID)
C POSY( J)= POSY( J)* LYREAL
CALL SSMUL( NUMCYL, LYREAL, POSY, XSTRID, POSY, XSTRID)
C IF( 0.5 .GE. DR( J)) THEN
C DR( J)= -1.0
C ELSE
C DR( J)= 1.0
C END IF
CALL SSGE( NUMCYL, 0.5, DR, XSTRID, LMASK1, YSTRID)
CALL SMASK( NUMCYL, -1.0, 0, 1.0, 0, LMASK1, XSTRID,
1 DR, YSTRID)
C DR( J)= DR( J)* ( A( J))**( 1.0/ ALPHA)
STEMP1= 1.0/ ALPHA
CALL SVPOW( NUMCYL, A, XSTRID, STEMP1, 0, VTEMP1, ZSTRID)
CALL SVMUL( NUMCYL, DR, XSTRID, VTEMP1, YSTRID, DR, XSTRID)
10 CONTINUE
C CARVE UP THE TOTAL INTERVAL INTO DL INCREMENTS.
DLX= LXREAL/ FLOAT( NUMXINC)
DLY= LYREAL/ FLOAT( NUMYINC)
C SET S. THIS CONTROLS THE DONUT DECAY.
S= 1
C PRIME THE MAX AND MIN.
IMGMAX= -1.0E30
IMGMIN= 1.0E30
C DIVIDE THE ROWS EVENLY ACROSS THE PROCESSORS.
DO 20 L= 1, NUMYINC, PRCSRS
LINCY= FLOAT( L- 1+ NODEID)* DLY
DO 25 J= 1, NUMXINC
LINCX= FLOAT( J- 1)* DLX
SUM= 0.0
C FOR PIXEL ( J, L), FIND THE DISTANCE TO EVERY CYLINDER CENTER:
C DO 30 N= 1, NUMCYL
C A( N)= SQRT(( LINCX- POSX( N))** 2+( LINCY- POSY( N))** 2)
CALL VPWLEDS( 3)
CALL SSSUB( NUMCYL, LINCX, POSX, XSTRID, VTEMP1, YSTRID)
CALL SVMUL( NUMCYL, VTEMP1, XSTRID, VTEMP1, XSTRID,
1 VTEMP1, XSTRID)
CALL SSSUB( NUMCYL, LINCY, POSY, XSTRID, VTEMP2, YSTRID)
CALL SVVTVP( NUMCYL, VTEMP2, XSTRID, VTEMP2, XSTRID,
1 VTEMP1, XSTRID, A, XSTRID)
CALL SVSQRT( NUMCYL, A, XSTRID, A, XSTRID)
C LMASK1( N)= A( N) .GE. RADIN( N) .AND. A( N) .LE. RADOUT( N)
CALL SGE( NUMCYL, A, XSTRID, RADIN, YSTRID, LMASK1, ZSTRID)
CALL SGT( NUMCYL, RADOUT, XSTRID, A, YSTRID, LMASK1, ZSTRID)
CALL LAND( NUMCYL, LMASK1, XSTRID, LMASK2, YSTRID, LMASK1, ZSTRID)
C IF( LMASK1( N)) THEN
C DEL AND SIGIN HAVE BEEN SQUARED AT GENERATION.
C SUM= SUM+ DR( N)*
C 1 EXP( - ABS( A( N)** 2- DEL( N))* SIGIN( N))**( 2* S)
C END IF
LENGTH= 0
CALL VPWLEDS( 0)
DO 130 N= 1, NUMCYL
C FIND THE DONUTS THAT MUST BE COMPUTED.
IF( LMASK1( N)) THEN
LENGTH= LENGTH+ 1
ITEMP1( LENGTH)= N
END IF
130 CONTINUE
CALL VPWLEDS( 1)
C GATHER THE VALUES OF THE VECTORS TO BE PROCESSED INTO CONTIGUOUS
C MEMORY.
IF( LENGTH .GT. 0) THEN
CALL SGATHR( LENGTH, A, XSTRID, ITEMP1, YSTRID,
1 A, XSTRID)
CALL SGATHR( LENGTH, DEL, XSTRID, ITEMP1, YSTRID,
1 VTEMP2, XSTRID)
CALL SGATHR( LENGTH, SIGIN, XSTRID, ITEMP1, YSTRID,
1 VTEMP3, XSTRID)
CALL SGATHR( LENGTH, DR, XSTRID, ITEMP1, YSTRID,
1 VTEMP1, XSTRID)
C NOTE: S IS ALWAYS 1.
C SUM= SUM+ DR( N)*
C 1 EXP( - ABS( A( N)** 2- DEL( N))* SIGIN( N))**( 2* S)
CALL VPWLEDS( 3)
CALL SVVTVM( LENGTH, A, XSTRID, A, XSTRID, VTEMP2, XSTRID,
1 A, XSTRID)
CALL SVABS( LENGTH, A, XSTRID, A, XSTRID)
CALL SSCAL( LENGTH, -1.0, A, XSTRID)
CALL SVMUL( LENGTH, A, XSTRID, VTEMP3, YSTRID, A, XSTRID)
CALL SVEXP( LENGTH, A, XSTRID, A, XSTRID)
CALL SVMUL( LENGTH, A, XSTRID, A, XSTRID, A, XSTRID)
SUM= SDOT( LENGTH, A, XSTRID, VTEMP1, YSTRID)
END IF
C30 CONTINUE
CR( J)= SUM
IF( SUM .GT. CRMAX) CRMAX= SUM
IF( SUM .LT. CRMIN) CRMIN= SUM
25 CONTINUE
C SEND THE LINE BACK TO THE HOST FOR STORAGE,
C ALONG WITH ROW NUMBER.
CR( NUMXINC+ 1)= FLOAT( L+ NODEID)
MSGTYP= 1
LENGTH= REAL4* CRBFSZ
CALL GREENLED( 1)
CALL VPWLEDS( 1)
CALL SENDW( CHANNL, MSGTYP, CR, LENGTH, HOST, HSTPID)
CALL VPWLEDS( 2)
CALL GREENLED( 0)
20 CONTINUE
CALL REDLED( 0)
C SEND THE MAX AND MIN FOR THIS NODE.
CR( 1)= CRMIN
CR( 2)= CRMAX
MSGTYP= 2
LENGTH= REAL4* 2
CALL SENDW( CHANNL, MSGTYP, CR, LENGTH, HOST, HSTPID)
CALL DETACH
STOP
END
SUBROUTINE FVRAND( A, C, N)
C MODULO RANDOM NUMBER GENERATOR FROM FPS.
INTEGER* 4 N
REAL* 4 A, C( N)
INTEGER* 4 I
REAL* 8 SEED, X
C EXTERNALS.
REAL* 8 DBLE, DMIN1
REAL* 4 SNGL
LOGICAL* 4 FIRST
DATA FIRST/ .TRUE./
IF( FIRST) THEN
SEED= DBLE( A)
FIRST= .FALSE.
END IF
DO 5 I= 1, N
X= SEED* 67108864.0D0
X= DMOD( 67081293.0D0* X+ 14181771.0D0, 67108864.0D0)
SEED= X/ 67108864.0
C( I)= SNGL( DMIN1( SEED, 0.999999D0))
5 CONTINUE
RETURN
END
Probably as much sense as the other crap put in this thread.
Has a nice Digital Look to me.