TABLE OF CONTENTS
- 1. NUM/save_1D
NUM/save_1D [ Modules ]
NOM
save_1D(calc)
DESCRIPTION
Posttraitement dans le cas ou le maillage est 1D
Post-processing if the mesh is 1D
ENTREES / INPUT
calc: objet de type calcul / object of type calcul
SORTIES / OUTPUT
SOURCE
! Ce logiciel est regi par la licence [CeCILL-B] ! This software is governed by the [CeCILL-B] license !=========================== DEBUT DES DECLARATIONS ==================== !.1----- Implicit, Use use num_typ IMPLICIT NONE !.2----- Declaration type(calcul), intent(inout) :: calc !! Objet calcul / Calcul object real(kind=kind(0.d0)) :: r,u,p,psiL,psiR,psi,xi,alf,err1,err2,volt character*7 :: fsol character*10 :: cn integer :: i,ic,id,iu,ios character(len=256) :: iomsg !=========================== DEBUT DU CODE EXECUTABLE ================== !------- Initialisations cn='0123456789' ic=int(calc%msh%numdom/100.) id=int((calc%msh%numdom-ic*100.)/10.) iu=calc%msh%numdom-id*10-ic*100 IF (calc%msh%numdom.gt.999) call print_err('save_1D','Domain index > 999 is not supported') fsol(5:5)=cn(ic+1:ic+1) fsol(6:6)=cn(id+1:id+1) fsol(7:7)=cn(iu+1:iu+1) ! ouverture des fichiers SELECT CASE (imodel_phy) CASE(2) fsol(1:4)='n_r_' ! densite OPEN(41,file=fsol,iostat=ios,iomsg=iomsg) if (ios /= 0) call print_err('save_1D','Cannot open '//trim(fsol)//' : '//trim(iomsg)) fsol(1:4)='n_u_' ! vitesse OPEN(42,file=fsol,iostat=ios,iomsg=iomsg) if (ios /= 0) call print_err('save_1D','Cannot open '//trim(fsol)//' : '//trim(iomsg)) fsol(1:4)='n_p_' ! pression OPEN(43,file=fsol,iostat=ios,iomsg=iomsg) if (ios /= 0) call print_err('save_1D','Cannot open '//trim(fsol)//' : '//trim(iomsg)) fsol(1:4)='n_e_' ! production numerique d entropie OPEN(49,file=fsol,iostat=ios,iomsg=iomsg) if (ios /= 0) call print_err('save_1D','Cannot open '//trim(fsol)//' : '//trim(iomsg)) fsol(1:4)='n_c_' ! fonction couleur OPEN(44,file=fsol,iostat=ios,iomsg=iomsg) if (ios /= 0) call print_err('save_1D','Cannot open '//trim(fsol)//' : '//trim(iomsg)) fsol(1:4)='e_r_' OPEN(51,file=fsol,iostat=ios,iomsg=iomsg) if (ios /= 0) call print_err('save_1D','Cannot open '//trim(fsol)//' : '//trim(iomsg)) fsol(1:4)='e_u_' OPEN(52,file=fsol,iostat=ios,iomsg=iomsg) if (ios /= 0) call print_err('save_1D','Cannot open '//trim(fsol)//' : '//trim(iomsg)) fsol(1:4)='e_p_' OPEN(53,file=fsol,iostat=ios,iomsg=iomsg) if (ios /= 0) call print_err('save_1D','Cannot open '//trim(fsol)//' : '//trim(iomsg)) fsol(1:4)='e_e_' open(54,file=fsol,iostat=ios,iomsg=iomsg) if (ios /= 0) call print_err('save_1D','Cannot open '//trim(fsol)//' : '//trim(iomsg)) CASE(1) fsol(1:4)='n_h_' ! hauteur d eau OPEN(41,file=fsol,iostat=ios,iomsg=iomsg) if (ios /= 0) call print_err('save_1D','Cannot open '//trim(fsol)//' : '//trim(iomsg)) fsol(1:4)='n_u_' ! vitesse OPEN(42,file=fsol,iostat=ios,iomsg=iomsg) if (ios /= 0) call print_err('save_1D','Cannot open '//trim(fsol)//' : '//trim(iomsg)) fsol(1:4)='n_e_' ! production numerique d entropie OPEN(49,file=fsol,iostat=ios,iomsg=iomsg) if (ios /= 0) call print_err('save_1D','Cannot open '//trim(fsol)//' : '//trim(iomsg)) fsol(1:4)='e_h_' OPEN(51,file=fsol,iostat=ios,iomsg=iomsg) if (ios /= 0) call print_err('save_1D','Cannot open '//trim(fsol)//' : '//trim(iomsg)) fsol(1:4)='e_u_' OPEN(52,file=fsol,iostat=ios,iomsg=iomsg) if (ios /= 0) call print_err('save_1D','Cannot open '//trim(fsol)//' : '//trim(iomsg)) open(53,file='cerf_1d.dat',iostat=ios,iomsg=iomsg) if (ios /= 0) call print_err('save_1D','Cannot open cerf_1d.dat : '//trim(iomsg)) write(53,*,iostat=ios,iomsg=iomsg)'TITLE="Solution Numerique CERF"' if (ios /= 0) call print_err('save_1D','Error writing cerf_1d header: '//trim(iomsg)) write(53,*,iostat=ios,iomsg=iomsg)' VARIABLES="X","h","U","zb","eta"' if (ios /= 0) call print_err('save_1D','Error writing cerf_1d variables: '//trim(iomsg)) write(53,111,iostat=ios,iomsg=iomsg)calc%msh%nb_cell if (ios /= 0) call print_err('save_1D','Error writing cerf_1d zone: '//trim(iomsg)) 111 format('ZONE T="", I=',i9,' F=POINT') CASE default call print_err('save_1D','Physical model not defined') END SELECT ! remplissage des fichiers de resultats num err1=0.d0 err2=0.d0 volt=0.d0 if (calc%tmax == 0.d0) call print_err('save_1D','tmax is zero, xi cannot be computed') DO i=1,calc%msh%nb_cell xi=calc%msh%list_cell(i)%center%x/calc%tmax CALL Bal2Prim(calc%msh%list_cell(i)%w) volt=volt+calc%msh%list_cell(i)%vol SELECT CASE (imodel_phy) CASE(2) write(41,*,iostat=ios,iomsg=iomsg) calc%msh%list_cell(i)%center%x,calc%msh%list_cell(i)%w%vprim%v(1)! densite if (ios /= 0) call print_err('save_1D','Error writing n_r_: '//trim(iomsg)) write(42,*,iostat=ios,iomsg=iomsg) calc%msh%list_cell(i)%center%x,calc%msh%list_cell(i)%w%vprim%v(2)! vitesse if (ios /= 0) call print_err('save_1D','Error writing n_u_: '//trim(iomsg)) write(43,*,iostat=ios,iomsg=iomsg) calc%msh%list_cell(i)%center%x,calc%msh%list_cell(i)%w%vprim%v(5)! pression if (ios /= 0) call print_err('save_1D','Error writing n_p_: '//trim(iomsg)) write(49,*,iostat=ios,iomsg=iomsg) calc%msh%list_cell(i)%center%x,calc%msh%list_cell(i)%w%vprim%v(6) if (ios /= 0) call print_err('save_1D','Error writing n_e_: '//trim(iomsg)) write(44,*,iostat=ios,iomsg=iomsg) calc%msh%list_cell(i)%center%x,calc%msh%list_cell(i)%w%vbal%v(5) if (ios /= 0) call print_err('save_1D','Error writing n_c_: '//trim(iomsg)) psiL=phy_state(1)%vbal%v(5) psiR=phy_state(2)%vbal%v(5) CALL riemisot(phy_state(1)%vprim%v(1),phy_state(1)%vprim%v(2),psiL,psiL, & phy_state(2)%vprim%v(1),phy_state(2)%vprim%v(2),psiR,psiR,xi,r,u,p,psi,alf) write(51,*,iostat=ios,iomsg=iomsg) calc%msh%list_cell(i)%center%x,r if (ios /= 0) call print_err('save_1D','Error writing e_r_: '//trim(iomsg)) write(52,*,iostat=ios,iomsg=iomsg) calc%msh%list_cell(i)%center%x,u if (ios /= 0) call print_err('save_1D','Error writing e_u_: '//trim(iomsg)) write(53,*,iostat=ios,iomsg=iomsg) calc%msh%list_cell(i)%center%x,p if (ios /= 0) call print_err('save_1D','Error writing e_p_: '//trim(iomsg)) write(54,*,iostat=ios,iomsg=iomsg) calc%msh%list_cell(i)%center%x,psi if (ios /= 0) call print_err('save_1D','Error writing e_e_: '//trim(iomsg)) CASE(1) write(41,*,iostat=ios,iomsg=iomsg) calc%msh%list_cell(i)%center%x,calc%msh%list_cell(i)%w%vprim%v(1)+& calc%msh%list_cell(i)%w%vprim%v(5)! hauteur d eau if (ios /= 0) call print_err('save_1D','Error writing n_h_: '//trim(iomsg)) write(42,*,iostat=ios,iomsg=iomsg) calc%msh%list_cell(i)%center%x,calc%msh%list_cell(i)%w%vprim%v(2)! vitesse if (ios /= 0) call print_err('save_1D','Error writing n_u_: '//trim(iomsg)) write(49,*,iostat=ios,iomsg=iomsg) calc%msh%list_cell(i)%center%x,calc%msh%list_cell(i)%w%vprim%v(4)! production entropie if (ios /= 0) call print_err('save_1D','Error writing n_e_: '//trim(iomsg)) CALL svriemann(phy_state(1)%vprim%v(1),phy_state(1)%vprim%v(2),& phy_state(2)%vprim%v(1),phy_state(2)%vprim%v(2),xi,r,u) write(51,*,iostat=ios,iomsg=iomsg) calc%msh%list_cell(i)%center%x,r if (ios /= 0) call print_err('save_1D','Error writing e_h_: '//trim(iomsg)) write(52,*,iostat=ios,iomsg=iomsg) calc%msh%list_cell(i)%center%x,u if (ios /= 0) call print_err('save_1D','Error writing e_u_: '//trim(iomsg)) write(53,*,iostat=ios,iomsg=iomsg) calc%msh%list_cell(i)%center%x,calc%msh%list_cell(i)%w%vprim%v(1),& calc%msh%list_cell(i)%w%vprim%v(2), & calc%msh%list_cell(i)%w%vprim%v(5), & calc%msh%list_cell(i)%w%vprim%v(1)+calc%msh%list_cell(i)%w%vprim%v(5) if (ios /= 0) call print_err('save_1D','Error writing cerf_1d.dat row: '//trim(iomsg)) CASE default call print_err('save_1D','Physical model not defined') END SELECT err1=err1+calc%msh%list_cell(i)%vol*(calc%msh%list_cell(i)%w%vprim%v(1)-r)**2 err2=err2+calc%msh%list_cell(i)%vol*(calc%msh%list_cell(i)%w%vprim%v(2)-u)**2 END DO if (volt <= 0.d0) call print_err('save_1D','Total cell volume is zero, cannot compute L2 errors') err1=sqrt(err1/volt) err2=sqrt(err2/volt) ! Fermeture des fichiers SELECT CASE (imodel_phy) CASE(2) write( *,1331)err1 write( *,1332)err2 1331 format('Density L2 error :',e13.5) 1332 format('Velocity L2 error :',e13.5) close(41) close(42) close(43) close(49) close(44) close(51) close(52) close(53) close(54) CASE(1) write( *,1001)err1 write( *,1002)err2 1001 format('Water level L2 error :',e13.5) 1002 format('Velocity L2 error :',e13.5) close(41) close(42) close(49) close(51) close(52) close(53) CASE default call print_err('save_1D','Unexpected physical model case') END SELECT !=========================== FIN DE LA ROUTINE ==================== END SUBROUTINE save_1D