Changeset 14841


Ignore:
Timestamp:
05/01/13 16:52:43 (12 years ago)
Author:
adhikari
Message:

CHG: some clean up

File:
1 edited

Legend:

Unmodified
Added
Removed
  • issm/trunk-jpl/src/c/modules/GiaDeflectionCorex/what0.f

    r14768 r14841  
    88     1decta(2),dyri1(nhank),dyri2(nhank),sna(nhank)
    99      double precision cinner(nhank),bcin(nhank)
    10 c     double precision trsl(Nrsl)
    1110      integer maxk
    12 c      double precision aswokm, asrpos, swok
    13 c      common /blocks/ aswokm,asrpos,swok
     11c ::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::
    1412      common /blockp/ pset
    1513      common /blockz/ zkp
    1614      common /blockm/ dekay1,dekay2,amp0,amp1,amp2
    17 c     common /blockk/ trsl
    1815      data yearco /3.15576d7/, pi /3.1415926535897932384d0/
    1916      data g /9.832186d0/, four /4.d0/, two /2.0d0/,
    2017     1 one /1.0d0/, zero/0.0d0/ , maxk/64/
     18c ::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::
    2119      twopi = two * pi
    2220      r2 = pset(6)
     
    2624      h  = pset(1)
    2725      urat = u1/u2
    28 c SA: ==================================================
     26c ::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::
     27c  alphap is dimensionless disk radius
     28c ::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::
    2929      alphap = pset(7)/pset(1)
    30 c     alphap is dimensionless disk radius
    3130      twoap = two * alphap
    3231      rghm = ( r1 * g * h * alphap ) / (two * u2)
     
    3837      dk = endk/dfloat(nhank)
    3938c
    40 c      if(icheck.eq.3) go to 999
    4139      if(idisk.gt.1) go to 7001
    4240      ak = zero
     
    5250      e2 = dexp(zkp2)
    5351      e4 = dexp(zkp4)
     52c ::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::
    5453      call freed(r2,u2,r1,u1,h,zkd,e1,e2,e4,b0,b1,a2,a1,a0,decay,amps)
    55 c
     54c ::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::
    5655      decta(1) = decay(1)/tmxyr
    5756      decta(2) = decay(2)/tmxyr
     
    5958      dyri2(ik) = decta(2)
    6059      sna(ik) = pikn
    61 c                            Form vectors for full construction in pwise.f
    62 c                            & stot.f
    63 c                            Note that freed will produce decay spectra
    64 c                            defined as positive, ie. negative decay must
    65 c                            reinsert a minus sign.
     60c ::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::
     61c Form vectors for full construction in pwise.f and stot.f
     62c Note that freed will produce decay spectra defined as positive,
     63c ie. negative decay must reinsert a minus sign.
     64c ::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::
    6665      dekay1(ik) = decay(1)
    6766      dekay2(ik) = decay(2)
     
    7372 7000 continue
    7473 7001 continue
    75 c     call dvecpr(zksam,nhank,'k sampled like Wolf',79,0,0)
    76 c     call dvecpr(zksamp,nhank,'k - dimensionless',79,0,0)
    77 c     call dvecpr(amp0,nhank,'amp0 of k from freed',79,0,0)
    78 c     call dvecpr(amp1,nhank,'amp1 of k from freed',79,0,0)
    79 c     call dvecpr(amp2,nhank,'amp2 of k from freed',79,0,0)
    80 c     call dvecpr(dekay1,nhank,'dekay1 of k from freed',79,0,0)
    81 c     call dvecpr(dekay2,nhank,'dekay2 of k from freed',79,0,0)
    82 c
    83 c A call below perfomed plots of the two decay times in years
    84 c sup {-1} vs. wave number k
    85 c
    86 c     if(icheck.gt.2) go to 49
    87 c     call wplot(dyri1,dyri2,sna)
    88 c     go to 999
    89 c
     74c ::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::
    9075c The following looped call sets up the free solution convolved with the
    9176c load function q hat.  Note that the returned vector set "cinner" is the
     
    9782c elliptical cross section.  Note loops 8500,8000 and 9500,9000 for the two
    9883c cases, respectively.
    99 c
     84c ::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::
    10085   49 go to (8499,9499), iedge
    101 c8499 tspan = trsl(irsl)
    102 c                                     trsl was nondimesionalized
    103 c                                     in the routine distme.f
    10486 8499 do 8000 ik = 1, nhank
    10587      xakap = zksamp(ik)*alphap
     
    11193      bcin(ik) = cinner(ik) * dbesj1(xakap)
    11294 8000 continue
    113 c      write(6,133) diku,xakap,pset(1),pset(7),alphap,urat
    114 c  133 format('      last diku  xakap in what0'/1h ,1p6d16.8)
    115 c      write(6,1) bcin(nhank),zksamp(nhank),fltng
    116 c    1 format('      last bcin  zksamp and fltng in what0'/1h ,1p3d16.8)
    117 c      write(6,2) cinner(nhank),cinner(1),qjadon
    118 c    2 format(' last and first cinner in what0 and qjadon'/1h ,1p3d16.8)
    119 c
     95c ::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::
    12096c "ojrule.f" computes the inverse Hankel trasform with a simple
    12197c Simpson's rule.  The routine "ojrule" is buliding a set of solutions stored
    12298c in common "blocks" in r or "asrpos(nrv) ", and computed rate or displacement
    12399c for each of N3G disks in "aswokm(nrv,N3G)" .
    124 c
     100c ::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::
    125101      call ojrule(dk,bcin)
    126102      go to 999
    127 c9499 tspan = trsl(irsl)
     103c ::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::
    128104 9499 do 9000 ik = 1, nhank
    129105      xakap = zksamp(ik)*alphap
     
    137113     1 - dcos(xakap) )
    138114 9000 continue
    139 c     write(6,1) bcin(nhank),zksamp(nhank)
    140115      call ojrule(dk,bcin)
    141 c     call dvecpr(cinner,nhank,'cinner at time t',79,0,0)
    142 c     call dvecpr(bcin,nhank,'bcin at time t',79,0,0)
    143116  999 return
     117
    144118      end
Note: See TracChangeset for help on using the changeset viewer.