Changeset 4050 for LMDZ6/trunk/libf/dyn3d/vlsplt.F
- Timestamp:
- Dec 23, 2021, 6:54:17 PM (3 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
LMDZ6/trunk/libf/dyn3d/vlsplt.F
r4008 r4050 4 4 5 5 SUBROUTINE vlsplt(q,pente_max,masse,w,pbaru,pbarv,pdt,iq) 6 USE infotrac, ONLY: nqtot, nqdesc,iqfils6 USE infotrac, ONLY: nqtot,tracers 7 7 c 8 8 c Auteurs: P.Le Van, F.Hourdin, F.Forget … … 83 83 CALL SCOPY(ijp1llm,masse,1,zm(1,1,iq),1) 84 84 85 if (nqdesc(iq).gt.0) then 86 do ifils=1,nqdesc(iq) 87 iq2=iqfils(ifils,iq) 88 CALL SCOPY(ijp1llm,q(1,1,iq2),1,zq(1,1,iq2),1) 89 enddo 90 endif !if (nqfils(iq).gt.0) then 85 do ifils=1,tracers(iq)%nqDescen 86 iq2=tracers(iq)%iqDescen(ifils) 87 CALL SCOPY(ijp1llm,q(1,1,iq2),1,zq(1,1,iq2),1) 88 enddo 91 89 92 90 cprint*,'Entree vlx1' … … 122 120 ENDDO 123 121 ! CRisi: aussi pour les fils 124 if (nqdesc(iq).gt.0) then 125 do ifils=1,nqdesc(iq) 126 iq2=iqfils(ifils,iq) 122 do ifils=1,tracers(iq)%nqDescen 123 iq2=tracers(iq)%iqDescen(ifils) 127 124 DO l=1,llm 128 DO ij=1,ip1jmp1129 q(ij,l,iq2)=zq(ij,l,iq2)130 ENDDO131 DO ij=1,ip1jm+1,iip1125 DO ij=1,ip1jmp1 126 q(ij,l,iq2)=zq(ij,l,iq2) 127 ENDDO 128 DO ij=1,ip1jm+1,iip1 132 129 q(ij+iim,l,iq2)=q(ij,l,iq2) 133 ENDDO130 ENDDO 134 131 ENDDO 135 enddo !do ifils=1,nqdesc(iq) 136 endif ! if (nqdesc(iq).gt.0) then 132 enddo 137 133 138 134 RETURN 139 135 END 140 136 RECURSIVE SUBROUTINE vlx(q,pente_max,masse,u_m,iq) 141 USE infotrac, ONLY : nqtot, nqfils,nqdesc,iqfils, ! CRisi137 USE infotrac, ONLY : nqtot,tracers, ! CRisi 142 138 & qperemin,masseqmin,ratiomin ! MVals et CRisi 143 139 … … 449 445 ! CRisi: appel récursif de l'advection sur les fils. 450 446 ! Il faut faire ça avant d'avoir mis à jour q et masse 451 !write(*,*) 'vlsplt 326: iq,nq fils(iq)=',iq,nqfils(iq)447 !write(*,*) 'vlsplt 326: iq,nqDesc(iq)=',iq,tracers(iq)%nqDescen 452 448 453 if (nqdesc(iq).gt.0) then 454 do ifils=1,nqdesc(iq) 455 iq2=iqfils(ifils,iq) 456 DO l=1,llm 449 do ifils=1,tracers(iq)%nqDescen 450 iq2=tracers(iq)%iqDescen(ifils) 451 DO l=1,llm 457 452 DO ij=iip2,ip1jm 458 ! On a besoin de q et masse seulement entre iip2 et ip1jm459 !masseq(ij,l,iq2)=masse(ij,l,iq)*q(ij,l,iq)460 !Ratio(ij,l,iq2)=q(ij,l,iq2)/q(ij,l,iq)461 !Mvals: veiller a ce qu'on n'ait pas de denominateur nul462 masseq(ij,l,iq2)=max(masse(ij,l,iq)*q(ij,l,iq),masseqmin)463 if (q(ij,l,iq).gt.qperemin) then464 Ratio(ij,l,iq2)=q(ij,l,iq2)/q(ij,l,iq)465 else466 Ratio(ij,l,iq2)=ratiomin467 endif453 ! On a besoin de q et masse seulement entre iip2 et ip1jm 454 !masseq(ij,l,iq2)=masse(ij,l,iq)*q(ij,l,iq) 455 !Ratio(ij,l,iq2)=q(ij,l,iq2)/q(ij,l,iq) 456 !Mvals: veiller a ce qu'on n'ait pas de denominateur nul 457 masseq(ij,l,iq2)=max(masse(ij,l,iq)*q(ij,l,iq),masseqmin) 458 if (q(ij,l,iq).gt.qperemin) then 459 Ratio(ij,l,iq2)=q(ij,l,iq2)/q(ij,l,iq) 460 else 461 Ratio(ij,l,iq2)=ratiomin 462 endif 468 463 enddo 469 enddo 470 enddo !do ifils=1,nqdesc(iq) 471 do ifils=1,nqfils(iq) 472 iq2=iqfils(ifils,iq) 473 call vlx(Ratio,pente_max,masseq,u_mq,iq2) 474 enddo !do ifils=1,nqfils(iq) 475 endif !if (nqfils(iq).gt.0) then 464 enddo 465 enddo 466 do ifils=1,tracers(iq)%nqDescen 467 iq2=tracers(iq)%iqDescen(ifils) 468 call vlx(Ratio,pente_max,masseq,u_mq,iq2) 469 enddo 476 470 ! end CRisi 477 471 … … 498 492 ! On calcule q entre iip2+1,ip1jm -> on fait pareil pour ratio 499 493 ! puis on boucle en longitude 500 if (nqfils(iq).gt.0) then 501 do ifils=1,nqdesc(iq) 502 iq2=iqfils(ifils,iq) 503 DO l=1,llm 494 do ifils=1,tracers(iq)%nqDescen 495 iq2=tracers(iq)%iqDescen(ifils) 496 DO l=1,llm 504 497 DO ij=iip2+1,ip1jm 505 498 q(ij,l,iq2)=q(ij,l,iq)*Ratio(ij,l,iq2) … … 508 501 q(ij-iim,l,iq2)=q(ij,l,iq2) 509 502 enddo ! DO ij=ijb+iip1-1,ije,iip1 510 enddo !DO l=1,llm 511 enddo !do ifils=1,nqdesc(iq) 512 endif !if (nqfils(iq).gt.0) then 503 enddo !DO l=1,llm 504 enddo 513 505 514 506 c CALL SCOPY((jjm-1)*llm,q(iip1+iip1,1),iip1,q(iip2,1),iip1) … … 519 511 END 520 512 RECURSIVE SUBROUTINE vly(q,pente_max,masse,masse_adv_v,iq) 521 USE infotrac, ONLY : nqtot, nqfils,nqdesc,iqfils, ! CRisi513 USE infotrac, ONLY : nqtot,tracers, ! CRisi 522 514 & qperemin,masseqmin,ratiomin ! MVals et CRisi 523 515 c … … 778 770 ! CRisi: appel récursif de l'advection sur les fils. 779 771 ! Il faut faire ça avant d'avoir mis à jour q et masse 780 !write(*,*) 'vly 689: iq,nq fils(iq)=',iq,nqfils(iq)772 !write(*,*) 'vly 689: iq,nqDesc(iq)=',iq,tracers(iq)%nqDescen 781 773 782 if (nqfils(iq).gt.0) then 783 do ifils=1,nqdesc(iq) 784 iq2=iqfils(ifils,iq) 785 DO l=1,llm 786 DO ij=1,ip1jmp1 787 ! attention, chaque fils doit avoir son masseq, sinon, le 1er 788 ! fils ecrase le masseq de ses freres. 789 !masseq(ij,l,iq2)=masse(ij,l,iq)*q(ij,l,iq) 790 !Ratio(ij,l,iq2)=q(ij,l,iq2)/q(ij,l,iq) 791 !MVals: veiller a ce qu'on n'ait pas de denominateur nul 792 masseq(ij,l,iq2)=max(masse(ij,l,iq)*q(ij,l,iq),masseqmin) 793 if (q(ij,l,iq).gt.qperemin) then 794 Ratio(ij,l,iq2)=q(ij,l,iq2)/q(ij,l,iq) 795 else 796 Ratio(ij,l,iq2)=ratiomin 797 endif 774 do ifils=1,tracers(iq)%nqDescen 775 iq2=tracers(iq)%iqDescen(ifils) 776 DO l=1,llm 777 DO ij=1,ip1jmp1 778 ! attention, chaque fils doit avoir son masseq, sinon, le 1er 779 ! fils ecrase le masseq de ses freres. 780 !masseq(ij,l,iq2)=masse(ij,l,iq)*q(ij,l,iq) 781 !Ratio(ij,l,iq2)=q(ij,l,iq2)/q(ij,l,iq) 782 !MVals: veiller a ce qu'on n'ait pas de denominateur nul 783 masseq(ij,l,iq2)=max(masse(ij,l,iq)*q(ij,l,iq),masseqmin) 784 if (q(ij,l,iq).gt.qperemin) then 785 Ratio(ij,l,iq2)=q(ij,l,iq2)/q(ij,l,iq) 786 else 787 Ratio(ij,l,iq2)=ratiomin 788 endif 798 789 enddo 799 enddo 800 enddo !do ifils=1,nqdesc(iq) 801 802 do ifils=1,nqfils(iq) 803 iq2=iqfils(ifils,iq) 804 call vly(Ratio,pente_max,masseq,qbyv,iq2) 805 enddo !do ifils=1,nqfils(iq) 806 endif !if (nqfils(iq).gt.0) then 790 enddo 791 enddo 792 793 do ifils=1,tracers(iq)%nqDescen 794 iq2=tracers(iq)%iqDescen(ifils) 795 call vly(Ratio,pente_max,masseq,qbyv,iq2) 796 enddo 807 797 808 798 DO l=1,llm … … 872 862 873 863 ! retablir les fils en rapport de melange par rapport a l'air: 874 if (nqfils(iq).gt.0) then 875 do ifils=1,nqdesc(iq) 876 iq2=iqfils(ifils,iq) 877 DO l=1,llm 864 do ifils=1,tracers(iq)%nqDescen 865 iq2=tracers(iq)%iqDescen(ifils) 866 DO l=1,llm 878 867 DO ij=1,ip1jmp1 879 868 q(ij,l,iq2)=q(ij,l,iq)*Ratio(ij,l,iq2) 880 869 enddo 881 enddo 882 enddo !do ifils=1,nqdesc(iq) 883 endif !if (nqfils(iq).gt.0) then 870 enddo 871 enddo 884 872 885 873 !write(*,*) 'vly 853: sortie' … … 888 876 END 889 877 RECURSIVE SUBROUTINE vlz(q,pente_max,masse,w,iq) 890 USE infotrac, ONLY : nqtot, nqfils,nqdesc,iqfils, ! CRisi878 USE infotrac, ONLY : nqtot,tracers, ! CRisi 891 879 & qperemin,masseqmin,ratiomin ! MVals et CRisi 892 880 c … … 1009 997 ! CRisi: appel récursif de l'advection sur les fils. 1010 998 ! Il faut faire ça avant d'avoir mis à jour q et masse 1011 !write(*,*) 'vlsplt 942: iq,nqfils(iq)=',iq,nqfils(iq) 1012 if (nqfils(iq).gt.0) then 1013 do ifils=1,nqdesc(iq) 1014 iq2=iqfils(ifils,iq) 1015 DO l=1,llm 999 !write(*,*) 'vlsplt 942: iq,nqChilds(iq)=',iq,nqChilds(iq) 1000 do ifils=1,tracers(iq)%nqDescen 1001 iq2=tracers(iq)%iqDescen(ifils) 1002 DO l=1,llm 1016 1003 DO ij=1,ip1jmp1 1017 !masseq(ij,l,iq2)=masse(ij,l,iq)*q(ij,l,iq)1018 !Ratio(ij,l,iq2)=q(ij,l,iq2)/q(ij,l,iq)1019 !MVals: veiller a ce qu'on n'ait pas de denominateur nul1020 masseq(ij,l,iq2)=max(masse(ij,l,iq)*q(ij,l,iq),masseqmin)1021 if (q(ij,l,iq).gt.qperemin) then1022 Ratio(ij,l,iq2)=q(ij,l,iq2)/q(ij,l,iq)1023 else1024 Ratio(ij,l,iq2)=ratiomin1025 endif1004 !masseq(ij,l,iq2)=masse(ij,l,iq)*q(ij,l,iq) 1005 !Ratio(ij,l,iq2)=q(ij,l,iq2)/q(ij,l,iq) 1006 !MVals: veiller a ce qu'on n'ait pas de denominateur nul 1007 masseq(ij,l,iq2)=max(masse(ij,l,iq)*q(ij,l,iq),masseqmin) 1008 if (q(ij,l,iq).gt.qperemin) then 1009 Ratio(ij,l,iq2)=q(ij,l,iq2)/q(ij,l,iq) 1010 else 1011 Ratio(ij,l,iq2)=ratiomin 1012 endif 1026 1013 enddo 1027 1028 enddo !do ifils=1,nqdesc(iq)1014 enddo 1015 enddo 1029 1016 1030 do ifils=1,nqfils(iq) 1031 iq2=iqfils(ifils,iq) 1032 call vlz(Ratio,pente_max,masseq,wq,iq2) 1033 enddo !do ifils=1,nqfils(iq) 1034 endif !if (nqfils(iq).gt.0) then 1017 do ifils=1,tracers(iq)%nqDescen 1018 iq2=tracers(iq)%iqDescen(ifils) 1019 call vlz(Ratio,pente_max,masseq,wq,iq2) 1020 enddo 1035 1021 ! end CRisi 1036 1022 … … 1045 1031 1046 1032 ! retablir les fils en rapport de melange par rapport a l'air: 1047 if (nqfils(iq).gt.0) then 1048 do ifils=1,nqdesc(iq) 1049 iq2=iqfils(ifils,iq) 1050 DO l=1,llm 1033 do ifils=1,tracers(iq)%nqDescen 1034 iq2=tracers(iq)%iqDescen(ifils) 1035 DO l=1,llm 1051 1036 DO ij=1,ip1jmp1 1052 1037 q(ij,l,iq2)=q(ij,l,iq)*Ratio(ij,l,iq2) 1053 1038 enddo 1054 enddo 1055 enddo !do ifils=1,nqdesc(iq) 1056 endif !if (nqfils(iq).gt.0) then 1039 enddo 1040 enddo 1057 1041 !write(*,*) 'vlsplt 1032' 1058 1042
Note: See TracChangeset
for help on using the changeset viewer.