C stecpu1:/home/g3/ogino/earthcg33/zbmain4902.f + zbsub4902.f C FILE NAME OGINO/EARTHB/OGNBB93 FILE NAME /POSC/OGNBB94 parameter (itapp=20+1*4,itappp=11) parameter (iip00=8,iiq00=4,iir00=150,thx=1.0) parameter(nnx=102,nny=102,nnz=2) c PARAMETER (X0=2.0,Y0=5.0,XL0=7.5,YL=7.5,ori2=360,last1=8) PARAMETER (X0=3.5,Y0=4.0,XL0=8.5,YL=8.5,ori2=360,last1=8) c PARAMETER (X0=1.2,Y0=5.0,XL0=5.5,YL=5.5) c PARAMETER (DXL=2.0,DYL=1.0,GTH=-45.0,TH0=60.0,DTH=10.0) PARAMETER (DXL=2.0,DYL=1.0,GTH=-45.0,TH0=60.0,DTH=10.0) c PARAMETER (NX= 320,NY= 80,NZ= 160) PARAMETER (NX= 180,NY= 60,NZ= 120) PARAMETER (N1=NX+2,N2=N1*(NY+2),N3=N2*(NZ+2)) c PARAMETER (NB=3,NBB=11,N4=N3*NB,N5=N3*NBB,THH0=70.0) PARAMETER (NB=3,NBB=11,N4=N3*NB,N5=N3*NBB,THH0=60.0) PARAMETER (NXG=40,NZG=20,ICU=1,MDFU=2) PARAMETER (MX=2*NX+3,JYG2=NZ,JYG=JYG2+2,MY=2*JYG2+3) PARAMETER (KK=MX*MY,MX2=MX,KP=0,LANK=20,KKK=N1*JYG*2,N6=KKK*4) c PARAMETER (LAST=1,IIQ0=1,XXL=160.50,YYL=40.50,ZZL=80.50,MM=2) PARAMETER (LAST= 1,IIQ0=1,XXL= 90.5,YYL=30.50,ZZL=60.50,MM=2) PARAMETER (NXP= 30,ARU=10.0,AR1=3.5,LAN1=30,LAN2=40,MOD=2) PARAMETER (IPEN=1,IAR=0,URMIN=0.01,BIS=-1.5E-4,EP1=1.0E-2) PARAMETER (RO01=5.0E-4,PR01=3.56E-8,VSW=0.044,VSWW=0.2*VSW) PARAMETER (JXG=N1*JYG,ARB=2.0,ARBI=3.5,N7=N3*3,NX1=NX+1) PARAMETER (NY1=NY+1,NZ1=NZ+1,NX2=NX+2,NY2=NY+2,NZ2=NZ+2) PARAMETER (N11=N1+1,N12=N2+1,N21=N1+N2,N22=N21+1,NXZ=NX*NZ) PARAMETER (NP1=NX2/2+NXP,NP1B=NX2/2-NXP) PARAMETER (NP2=NP1*NY2,NP3=NP2*NZ2,NP4=NP3*3) PARAMETER (IBB=30,IBL=120,IBD=30,KBB=2,KBL=14,KBD=4) PARAMETER (M0=1,M1=N2/M0,M2=NZ2*8*M0,M3=NZ2*5*M0) parameter (idci=2,ity=2,isz=1,fct=1.00,H1=0.45,H2=0.35) c c parameter (ori1= 45.0,tim1= 720.0,XLIM=0.0,XLI2=10.0) c parameter (ori1= 50.0,tim1= 960.0,XLIM=0.0,XLI2=10.0) c parameter (ori1= 55.0,tim1=1200.0,XLIM=0.0,XLI2=10.0) c parameter (ori1= 60.0,tim1= 960.0,XLIM=0.0,XLI2=10.0) c parameter (ori1= 65.0,tim1=1200.0,XLIM=0.0,XLI2=10.0) c parameter (ori1= 70.0,tim1=1200.0,XLIM=0.0,XLI2=10.0) c parameter (ori1= 75.0,tim1= 960.0,XLIM=0.0,XLI2=10.0) parameter (ori1= 80.0,tim1=1200.0,XLIM=0.0,XLI2=10.0) c parameter (ori1= 85.0,tim1= 960.0,XLIM=0.0,XLI2=10.0) c parameter (ori1= 90.0,tim1= 720.0,XLIM=0.0,XLI2=10.0) c c character tl1*40/'Northward turning from southward IMF'/ c character tl2*40/'Bz=0nT -5nT 5nT t=161.5m (104.5m)'/ c character tl1*40/'Southward turning from northward IMF'/ c character tl2*40/'Bz=0nT 5nT -5nT t=171.0m (114.0m)'/ c c character chrt*32/"Rotation of Incoming IMF "/ character chrt*32/"Event for November 17, 1996 "/ c character chrt*32/"Event for March 19, 1999 "/ c character chrt*32/"Incoming uniform IMF "/ c character chrt*32/"No uniform IMF "/ character chg(3)*46 c data chg(1)/"MHD Simulation for 1996 November 17 Event "/ c data chg(1)/"MHD Simulation for 1999 March 19 Event "/ c data chg(1)/"MHD Simulation for 1999 October 21 Event "/ data chg(1)/"3D MHD Simulation of Earth's Magnetosphere "/ c data chg(1)/"MHD Simulation for 1999 October 22 Event "/ c data chg(2)/"B= 2.5nT t=960m (240m) "/ data chg(3)/"Density and energy of cross section "/ c character chr(10)*46 data chr(1)/"T = 21-10-1999 00:11:00 "/ data chr(2)/"T = 21-10-1999 00:12:00 "/ data chr(3)/"T = 21-10-1999 00:13:00 "/ data chr(4)/"T = 21-10-1999 00:14:00 "/ data chr(5)/"T = 21-10-1999 00:15:00 "/ data chr(6)/"T = 21-10-1999 00:16:00 "/ data chr(7)/"T = 21-10-1999 00:17:00 "/ data chr(8)/"T = 21-10-1999 00:18:00 "/ data chr(9)/"T = 21-10-1999 00:19:00 "/ data chr(10)/"T = 21-10-1999 00:20:00 "/ c DIMENSION F(N4),P(N6) DIMENSION IZM(9),ZMIN(9),ZMAX(9),IA(10),AA(24) common app(nnx,nny,nnz),ie(nnx,nny,nnz) DATA IZM/2,2,2,2,2,2,2,2,2/ DATA ZMIN/-3.0E-2, 0.0,-0.0E-3, 0.0, -0.02, & -0.5E-3,-0.4E-4,-0.5E-3, 0.0/ DATA ZMAX/ 3.0E-2, 0.5E-6, 2.0E-3, 0.5E-3, 0.02, & 0.5E-3, 0.4E-4, 0.5E-3, 0.5E-5/ c c open(11,file='eartha10.data', 1 access='sequential',form='unformatted') c c call xyopen ino=0 DO 300 JJJ=1,1 ITAP=JJJ+10 IA(1)=NX IA(2)=JYG2 IA(3)=NX IA(4)=NY IA(5)=NZ IA(6)=NXP IA(7)=6 IA(8)=2 AA(1)=TH0 AA(2)=ARU AA(3)=AR1 AA(4)=ARB AA(5)=XXL AA(6)=YYL AA(7)=ZZL AA(8)=1.0 AA(9)=0.2 AA(11)=-1.0 AA(12)=0.0 AA(15)=XL0 AA(16)=YL AA(17)=GTH AA(18)=-20.0 c AA(19)=120.0 AA(19)=60.0 AA(21)=25.0 AA(20)=EP1 AA(23)=XLIM AA(24)=XLI2 XEP=0.5*FLOAT(NX1-2*NXP)/FLOAT(NX1) C C GRAPHIC OPEN c CALL XYOPEN(0,0.0,0.0,33.0,33.0,0) XLL1=3.0*XL+2.0*DXL+X0 AR2=AR1*AR1 PI=3.1415926 NB1=NB+1 NXP2=NXP*2 MXA=MX-2 MYA=MY-2 HX=XXL/FLOAT(NX1) HY=YYL/FLOAT(NY1) HZ=ZZL/FLOAT(NZ1) HX2=0.5*HX HY2=0.5*HY HZ2=0.5*HZ HXG=FLOAT(NX1)/FLOAT(NXG+1) HZG=FLOAT(JYG2+1)/FLOAT(NZG+1) HXX=XXL/FLOAT(NXG+1) HZZ=(YYL+ZZL-HY2-HZ2)/FLOAT(NZG+1) CC HZZ=(ZZL-HY2-HZ2)/FLOAT(NZG+1) NXZG=NXG*NZG C IA(1)=NXG IA(2)=NZG AA(1)=THH0 AA(4)=ARBI AD=2.0 TH1=GTH*PI/180.0 XC1=COS(TH1) XS1=SIN(TH1) YLL1=YL+DYL XL=XL0*0.5*(AA(19)-AA(18))/AA(21) AA(15)=XL AA(16)=YL xla=xl*(-aa(18))/(aa(19)-aa(18)) yla=yl*0.5 xlb=yla*xc1 ylb=-yla*xs1 IIQ=IIQ0-1 JJ=0 C DO 1003 II=1,LAST IIQ=IIQ+1 DO 724 I1=1,M2 I2=M1*(I1-M3-1)+1 I3=M1*(I1-M3) IF(I1.LE.M3) READ(ITAP) (F(I),I=1,M1) IF(I1.GT.M3) READ(ITAP) (F(I),I=I2,I3) 724 CONTINUE c c ll=3 do 172 m=1,ll do 172 j=1,ny2 c do 174 k=1,nz2 do 174 i=1,nx2 i1=i+n1*(j-1)+n2*(k-1)+n3*(m-1) j1=i+n1*(nz2-k) p(j1)=f(i1) 174 continue c do 172 k=1,nz2 do 172 i=1,nx2 i1=i+n1*(j-1)+n2*(k-1)+n3*(m-1) j1=i+n1*(k-1) f(i1)=p(j1) if(m.eq.1) f(i1)=-f(i1) if(m.eq.2) f(i1)=-f(i1) c if(m.eq.7) f(i1)=-f(i1) 172 continue c c cc do 555 i=m1*2*nz2*2+1,m1*2*nz2*2.5 cc F(i)=-F(i) cc 555 continue IF(IIQ.NE.IIQ0) GO TO 1003 IIQ=0 IF(II.LE.0) GO TO 1003 C X1=2.0 C CALL BOUND3(NX,NY,NZ,NXP,MDFU,NB,BIS,RO01,PR01,VSW, C & X1,HX,HY,HZ,ARU,AR2,F,P) C c CALL PLOTS(NAME,16.2) c CALL FACTOR(1.50) c CALL PLOT(0.0,0.0,-3) c ino=ino+1 call plots(idci,ity,isz,ino) call factor(fct) CALL NEWPEN(IPEN) I1=II c CALL DATA(0.5,0.5,LAST,I1,NXP,NX) C C 3 B BAGNETIC FIELD BX,BY,BZ IA(1)=NXG/2 IA(2)=NZG/2 IA(1)=NXG IA(2)=NXG c IA(1)=100 c IA(2)=100 IA(1)=50 IA(2)=50 IM=3 KP1=IZM(IM) VMIN=ZMIN(IM) VMAX=ZMAX(IM) X1=X0-0.5*YL*XC1 AA(13)=AMAX1(X0,X1) c AA(14)=Y0-0.05*YL AA(14)=Y0+0.00*YL+0.00*DYL IA(7)=6 write(6,*) ino,' No. 01' cc polar cap start B Vector under figure c CALL AINTE1A(IA,AA,F,P) c write(6,*) ino,' No. 012' call zsub33(ia,aa) c write(6,*) ino,' No. 013' c CALL AINTE1(IA,AA,F,P) xb=x0 yb=aa(14) call symbol(xb-h2*2.0,yb+yla-h2*0.5,h2,'X',0.0,1) call symbol(xb+xla+xlb-h2*0.5,yb+yla-ylb-h2*1.5,h2,'Y',0.0,1) call symbol(xb+xla-h2*0.5,yb+yl+h2*0.5,h2,'Z',0.0,1) C IA(1)=NXG IA(1)=NXG/2 IA(2)=NZG IA(2)=10 AA(13)=AMAX1(X0,X1) c AA(14)=Y0+1.30*YL+DYL AA(14)=Y0+0.00*YL+0.00*DYL IA(7)=6 write(6,*) ino,' No. 02' cc from equator start B Vector upper figure c CALL AINTE2(IA,AA,F,P) CALL AINTE21(IA,AA,F,P) yb=aa(14) call symbol(xb-h2*2.0,yb+yla-h2*0.5,h2,'X',0.0,1) call symbol(xb+xla+xlb-h2*0.5,yb+yla-ylb-h2*1.5,h2,'Y',0.0,1) call symbol(xb+xla-h2*0.5,yb+yl+h2*0.5,h2,'Z',0.0,1) call symbol(xb-h2*4.0,yb+yla-h2*2.0,h2,'20Re',0.0,4) c call symbol(xb+xl+h2*3.0,yb+yla-h2*2.0,h2,'-120Re',0.0,6) call symbol(xb+xl+h2*3.0,yb+yla-h2*2.0,h2,'- 60Re ',0.0,6) call symbol(xb+xla+xlb+h2*2.0,yb+yla-ylb-h2*1.5,h2,'25Re',0.0,4) call symbol(xb+xla+h2*2.0,yb+yl+h2*0.5,h2,'25Re',0.0,4) c xgb=x0+yl*0.2 xgg=xgb+h1*20.0 ygg=y0+yl*1.10+dyl*0.5+h1*0.5 ccccccccccccccccccccccccccccc ori=ori1 c 500 call kakudo(xgg,ygg+h1*3.5,ori,h1*1.0) c call symblc(xgg+5.7*h1,ygg+h1*3.5,h1*1.0,"\\260",0.0,5) c call symblc(xgg-h1*1.3,ygg+h1*3.5,h1*1.0,"\\161",0.0,5) c hh2=h1*6 c tim=tim1 c xg1=xgb+hh2+h1*2 xg1=xgb+hh2+h1*4 yg1=ygg+h1*3.5 c call timesg(xg1,yg1,tim,h1) c call symbol(xg1-h1*2,yg1,h1*1.0," t=",0.0,3) c call symbol(xg1+h1*5,yg1,h1*1.0,"min",0.0,3) c call symbol(xgb,yg1,h1*1.0,chg,0.0,8) call symbol(xgb,yg1,h1*1.0,chr(jjj),0.0,42) call symbol(xgb,yg1+h1*2.5,h1*1.2,chg(1),0.0,46) c call symbol(xb,yg1,h1,tl2,0.0,40) c call symblc(xb+h1*5.4,yg1+h1*0.2,h1*0.7,'\256',0.0,4) c call symblc(xb+h1*10.2,yg1+h1*0.2,h1*0.7,'\256',0.0,4) c call symblc(xb+h1*09.7,yg1+h1*0.2,h1*0.7,'\256',0.0,4) c call symbol(xb,yg1+h1*1.0,h1,tl1,0.0,40) CALL PLOTE 1003 CONTINUE 300 CONTINUE C 999 CONTINUE c CALL XYCLOS STOP END subroutine kakudo(xg,yg,ori,h1) character isymb*20 write(isymb,40) ori 40 format("=",f6.1) call symbol(xg,yg,h1,isymb,0.0,8) c write(6,*) "(",isymb,") s" return end c subroutine timesg(xg,yg,ori,h1) character isymb*20 write(isymb,40) ori 40 format(" ",f6.1) call symbol(xg,yg,h1,isymb,0.0,7) c write(6,*) "(",isymb,") s" return end