source: PSPA/parmelaPSPA/trunk/trwave.f @ 412

Last change on this file since 412 was 12, checked in by lemeur, 12 years ago

parmela pspa initial

File size: 3.2 KB
Line 
1       subroutine trwave(zz,rr,wez,wer,wbphi,i,wtp)
2c       nn is used in this subroutine, so common/misc/
3c       cannot be used here.
4c---get field of traveling wave at r,z, and time wtp
5c----------------------------------------------------------------------
6c
7      include 'param_sz.h'
8      include 'cfldscom.h'
9      include 'constcom.h'
10      include 'coordcom.h'
11      include 'debug.h'
12      include 'pcordcom.h'
13      include 'syscom.h'
14      include 'tstepcom.h'
15      include 'wavescom.h'
16      include 'ucom.h'
17c       
18c--------------------------------------------------------------------------
19c*
20      wer=0.
21      wez=0.
22      wbphi=0.
23c------------------------- if nprint gt 0 ---------------------------------
24      if (rne.eq.0.)then
25        do 100 i=1,nel
26        if(zz.ge.zloc(i-1).and.zz.le.zloc(i).and.ntype(i).eq.9)goto 200
27c---if zz is in range of a cell card output the field of that cell
28        if(zz.ge.zloc(i-1).and.zz.le.zloc(i).and.ntype(i).eq.7)then
29           ne=i
30           rne=ne
31           ph=wtp*radian*el(4,ne)/freq
32           sp=sin(ph)
33           cp=cos(ph)
34           nc=el(6,ne)
35           zc=zz+hcll(nc)-zloc(ne)
36           if (el(8,ne).gt.0.) zc=zc-hcll(nc)
37           call cfield(rr,zc,nc,ez,er,bt)
38           s=sp*el(10,ne)+cp*el(9,ne)
39           wez=el(5,ne)*100.*ez*s
40           c=cp*el(10,ne)-sp*el(9,ne)
41           wer=el(5,ne)*100.*er*s
42           wbphi=el(5,ne)*100.*bt*c
43           rne=0.
44        endif ! cell
45 100    continue
46        i=0
47        return
48c---the z position is within the defined region of the i-th traveling wave---
49 200    continue
50        ne=i
51        nc=el(6,ne)
52        i=ne+1-nelw(nc)
53c--------------------- end of if nprint gt 0 --------------------------------
54      else
55        ne=rne
56        nc=el(6,ne)
57        i=ne+1-nelw(nc)
58      endif
59c---------now add up all the spatial components from nmin to nmax-----------
60      amp0=el(5,ne)
61      do 300 n=nmin,nmax
62        nn=n-nmin+1
63        cay=twopi/wavel*el(8,ne)/freq
64        cay3=wk3(nn,i,nc)
65        cay1=wk1(nn,i,nc)
66cliu    cay1rsq=(cay1*rr)**2
67        cay1rsq=-cay1*rr**2         
68        if(ne.ne.nelw(nc))then
69          phase1=wtp*radian*el(8,ne)/freq
70c       ------------------------------------------
71          if(ne.ne.nw(nc)) then
72                phase2=-cay3*(zz-zloc(ne-1)-el(1,ne)/2.)
73          else
74                phase2=-cay3*(zz-zloc(ne-1)-el(1,ne))
75          endif
76c       ----------------------------------------------
77          phase3=wpoff(nn,i,nc)
78          phase4=el(4,ne)
79          phase=phase1+phase2+phase3+phase4
80c         phase=wtp*radian*el(8,ne)/freq-cay3*(zz-zloc(ne-1)
81C     &         -el(1,ne)/2.)+wpoff(nn,i,nc)+el(4,ne)
82        else
83          phase1=wtp*radian*el(8,ne)/freq
84          phase2=-cay3*(zz-zloc(ne-1))
85          phase3=wpoff(nn,i,nc)
86          phase4=el(4,ne)
87          phase=phase1+phase2+phase3+phase4
88c         phase=wtp*radian*el(8,ne)/freq-cay3*(ZZ-zloc(ne-1))
89C     &         +wpoff(nn,i,nc)+el(4,ne)
90        endif
91         amplitd=amp0*wa(nn,i,nc)
92         wez= wez+amplitd*xi0(cay1rsq)*cos(phase)
93         dern= -amplitd*xi1ox(cay1rsq)*rr
94cliu     wer=wer + dern*cos(phase)*cay3
95         wer=wer + dern*sin(phase)*cay3   
96         wbphi= wbphi+ dern*sin(phase)*cay
97 300  continue
98      return
99      end
100c++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++*
Note: See TracBrowser for help on using the repository browser.