source: trunk/LMDZ.GENERIC/libf/phystd/stellarlong.F @ 1351

Last change on this file since 1351 was 1308, checked in by emillour, 10 years ago

Generic GCM:
Some cleanup to simplify dynamics/physics interactions by getting rid
of dimphys.h (i.e. the nlayermx parameter) and minimizing use of
dimension.h in the physics.
EM

File size: 2.5 KB
Line 
1      SUBROUTINE stellarlong(pday,pstellong)
2     
3      USE planete_mod, ONLY: year_day, peri_day, e_elips, timeperi
4      IMPLICIT NONE
5
6c=======================================================================
7c
8c   Objet:
9c   ------
10c
11c      Calcul de la distance soleil-planete et de la declinaison
12c   en fonction du jour de l'annee.
13c
14c
15c   Methode:
16c   --------
17c
18c      Calcul complet de l'elipse
19c
20c   Interface:
21c   ----------
22c
23c      Uncommon comprenant les parametres orbitaux.
24c
25c   Arguments:
26c   ----------
27c
28c   Input:
29c   ------
30c   pday          jour de l'annee (le jour 0 correspondant a l'equinoxe)
31c
32c   Output:
33c   -------
34c   pdist_star     distance entre le soleil et la planete
35c                 ( en unite astronomique pour utiliser la constante
36c                  solaire terrestre 1370 Wm-2 )
37c   pdecli        declinaison ( en radians )
38c
39c=======================================================================
40c-----------------------------------------------------------------------
41c   Declarations:
42c   -------------
43
44!#include "planete.h"
45#include "comcstfi.h"
46
47c arguments:
48c ----------
49
50      REAL pday,pdist_star,pdecli,pstellong
51      LOGICAL lwrite
52
53c Local:
54c ------
55
56      REAL zanom,xref,zx0,zdx,zteta,zz
57      INTEGER iter
58
59
60c-----------------------------------------------------------------------
61c calcul de l'angle polaire et de la distance au soleil :
62c -------------------------------------------------------
63
64c  calcul de l'zanomalie moyenne
65
66      zz=(pday-peri_day)/year_day
67      pi=2.*asin(1.)
68      zanom=2.*pi*(zz-nint(zz))
69      xref=abs(zanom)
70
71c  resolution de l'equation horaire  zx0 - e * sin (zx0) = xref
72c  methode de Newton
73
74      zx0=xref+e_elips*sin(xref)
75      DO 110 iter=1,10
76         zdx=-(zx0-e_elips*sin(zx0)-xref)/(1.-e_elips*cos(zx0))
77         if(abs(zdx).le.(1.e-7)) goto 120
78         zx0=zx0+zdx
79110   continue
80120   continue
81      zx0=zx0+zdx
82      if(zanom.lt.0.) zx0=-zx0
83
84c zteta est la longitude solaire
85
86      zteta=2.*atan(sqrt((1.+e_elips)/(1.-e_elips))*tan(zx0/2.))
87
88      pstellong=zteta-timeperi
89
90      IF(pstellong.LT.0.) pstellong=pstellong+2.*pi
91      IF(pstellong.GT.2.*pi) pstellong=pstellong-2.*pi
92c-----------------------------------------------------------------------
93c   sorties eventuelles:
94c   ---------------------
95
96c     IF (lwrite) THEN
97c        PRINT*,'jour de l"annee   :',pday
98c        PRINT*,'distance au soleil (en unite astronomique) :',pdist_star
99c        PRINT*,'declinaison (en degres) :',pdecli*180./pi
100c     ENDIF
101
102      RETURN
103      END
Note: See TracBrowser for help on using the repository browser.