Skip to content

Commit ecaf888

Browse files
cpintecpinte
committed
Simplifying dust_prop
1 parent 3c7888e commit ecaf888

File tree

1 file changed

+69
-188
lines changed

1 file changed

+69
-188
lines changed

‎src/dust_prop.f90‎

Lines changed: 69 additions & 188 deletions
Original file line numberDiff line numberDiff line change
@@ -817,7 +817,7 @@ subroutine opacite(lambda, p_lambda, no_scatt)
817817

818818
! Attention : dans le cas no_strat, il ne faut pas que la cellule (1,1,1) soit vide.
819819
! on la met à nbre_grains et on effacera apres
820-
! c'est pour les prop de diffusion en relatif donc la veleur exacte n'a pas d'importante
820+
! c'est pour les prop de diffusion en relatif donc la valeur exacte n'a pas d'importante
821821
ldens0 = .false.
822822
if (.not.lvariable_dust) then
823823
if (icell_not_empty > icell1) then ! icell1==1
@@ -841,80 +841,38 @@ subroutine opacite(lambda, p_lambda, no_scatt)
841841
kappa(icell,lambda) = 0.0
842842
k_sca_tot = 0.0
843843

844-
if (lvariable_dust) then
845-
do k=1,n_grains_tot ! Expensive when n_cells is large
846-
density = dust_density(k,icell) * nbre_grains(k)
847-
kappa(icell,lambda) = kappa(icell,lambda) + C_ext(k,lambda) * density
848-
k_sca_tot = k_sca_tot + C_sca(k,lambda) * density
849-
enddo !k
850-
else
851-
d1 = dust_density(1,icell)
852-
if (d1 > tiny_dp) then
853-
do k=1,n_grains_tot ! Expensive when n_cells is large
854-
density = d1 * nbre_grains(k)
855-
kappa(icell,lambda) = kappa(icell,lambda) + C_ext(k,lambda) * density
856-
k_sca_tot = k_sca_tot + C_sca(k,lambda) * density
857-
enddo !k
858-
endif
859-
endif
844+
do k=1,n_grains_tot ! Expensive when n_cells is large
845+
density = dust_density(k,icell) * nbre_grains(k)
846+
kappa(icell,lambda) = kappa(icell,lambda) + C_ext(k,lambda) * density
847+
k_sca_tot = k_sca_tot + C_sca(k,lambda) * density
848+
enddo !k
860849

861850
if (kappa(icell,lambda) > tiny_real) tab_albedo_pos(icell,lambda) = k_sca_tot/kappa(icell,lambda)
862851

863852
if (aniso_method==2) then
864853
tab_g_pos(icell,lambda) = 0.0
865-
if (lvariable_dust) then
866-
do k=1,n_grains_tot ! Expensive when n_cells is large
867-
density=dust_density(k,icell) * nbre_grains(k)
868-
tab_g_pos(icell,lambda) = tab_g_pos(icell,lambda) + &
869-
C_sca(k,lambda) * density * tab_g(k,lambda)
870-
enddo ! k
871-
else
872-
d1 = dust_density(1,icell)
873-
if (d1 > tiny_dp) then
874-
do k=1,n_grains_tot ! Expensive when n_cells is large
875-
density=d1 * nbre_grains(k)
876-
tab_g_pos(icell,lambda) = tab_g_pos(icell,lambda) + &
877-
C_sca(k,lambda) * density * tab_g(k,lambda)
878-
enddo ! k
879-
endif
880-
endif
854+
do k=1,n_grains_tot ! Expensive when n_cells is large
855+
density=dust_density(k,icell) * nbre_grains(k)
856+
tab_g_pos(icell,lambda) = tab_g_pos(icell,lambda) + &
857+
C_sca(k,lambda) * density * tab_g(k,lambda)
858+
enddo ! k
881859
if (k_sca_tot > tiny_real) tab_g_pos(icell,lambda) = tab_g_pos(icell,lambda)/k_sca_tot
882860
endif
883861

884862
if (lRE_LTE) then
885863
kappa_abs_LTE(icell,lambda) = 0.0
886-
if (lvariable_dust) then
887-
do k=grain_RE_LTE_start,grain_RE_LTE_end ! Expensive when n_cells is large
888-
kappa_abs_LTE(icell,lambda) = kappa_abs_LTE(icell,lambda) + &
889-
C_abs(k,lambda) * dust_density(k,icell) * nbre_grains(k)
890-
enddo
891-
else
892-
d1 = dust_density(1,icell)
893-
if (d1 > tiny_dp) then
894-
do k=grain_RE_LTE_start,grain_RE_LTE_end ! Expensive when n_cells is large
895-
kappa_abs_LTE(icell,lambda) = kappa_abs_LTE(icell,lambda) + &
896-
C_abs(k,lambda) * d1 * nbre_grains(k)
897-
enddo
898-
endif
899-
endif
864+
do k=grain_RE_LTE_start,grain_RE_LTE_end ! Expensive when n_cells is large
865+
kappa_abs_LTE(icell,lambda) = kappa_abs_LTE(icell,lambda) + &
866+
C_abs(k,lambda) * dust_density(k,icell) * nbre_grains(k)
867+
enddo
900868
endif
901869

902870
if (lRE_nLTE) then
903871
kappa_abs_nLTE(icell,lambda) = 0.0
904-
if (lvariable_dust) then
905-
do k=grain_RE_nLTE_start,grain_RE_nLTE_end
906-
kappa_abs_nLTE(icell,lambda) = kappa_abs_nLTE(icell,lambda) + &
907-
C_abs(k,lambda) * dust_density(k,icell) * nbre_grains(k)
908-
enddo
909-
else
910-
d1 = dust_density(1,icell)
911-
if (d1 > tiny_dp) then
912-
do k=grain_RE_nLTE_start,grain_RE_nLTE_end
913-
kappa_abs_nLTE(icell,lambda) = kappa_abs_nLTE(icell,lambda) + &
914-
C_abs(k,lambda) * d1 * nbre_grains(k)
915-
enddo
916-
endif
917-
endif
872+
do k=grain_RE_nLTE_start,grain_RE_nLTE_end
873+
kappa_abs_nLTE(icell,lambda) = kappa_abs_nLTE(icell,lambda) + &
874+
C_abs(k,lambda) * dust_density(k,icell) * nbre_grains(k)
875+
enddo
918876
endif
919877
enddo ! icell
920878

@@ -926,32 +884,17 @@ subroutine opacite(lambda, p_lambda, no_scatt)
926884
k_abs_tot = 0.0
927885
k_abs_RE = 0.0
928886
k_abs_LTE = 0.0
929-
if (lvariable_dust) then
930-
do k=1, n_grains_tot
931-
k_abs_tot = k_abs_tot + C_abs(k,lambda) * dust_density(k,icell) * nbre_grains(k)
932-
enddo
933-
do k=grain_RE_LTE_start,grain_RE_LTE_end
934-
k_abs_LTE = k_abs_LTE + C_abs(k,lambda) * dust_density(k,icell) * nbre_grains(k)
935-
enddo
936-
k_abs_RE = k_abs_LTE
937-
do k=grain_RE_nLTE_start,grain_RE_nLTE_end
938-
k_abs_RE = k_abs_RE + C_abs(k,lambda) * dust_density(k,icell) * nbre_grains(k)
939-
enddo
940-
else
941-
d1 = dust_density(1,icell)
942-
if (d1 > tiny_dp) then
943-
do k=1, n_grains_tot
944-
k_abs_tot = k_abs_tot + C_abs(k,lambda) * d1 * nbre_grains(k)
945-
enddo
946-
do k=grain_RE_LTE_start,grain_RE_LTE_end
947-
k_abs_LTE = k_abs_LTE + C_abs(k,lambda) * d1 * nbre_grains(k)
948-
enddo
949-
k_abs_RE = k_abs_LTE
950-
do k=grain_RE_nLTE_start,grain_RE_nLTE_end
951-
k_abs_RE = k_abs_RE + C_abs(k,lambda) * d1 * nbre_grains(k)
952-
enddo
953-
endif
954-
endif
887+
888+
do k=1, n_grains_tot
889+
k_abs_tot = k_abs_tot + C_abs(k,lambda) * dust_density(p_k,icell) * nbre_grains(k)
890+
enddo
891+
do k=grain_RE_LTE_start,grain_RE_LTE_end
892+
k_abs_LTE = k_abs_LTE + C_abs(k,lambda) * dust_density(p_k,icell) * nbre_grains(k)
893+
enddo
894+
k_abs_RE = k_abs_LTE
895+
do k=grain_RE_nLTE_start,grain_RE_nLTE_end
896+
k_abs_RE = k_abs_RE + C_abs(k,lambda) * dust_density(p_k,icell) * nbre_grains(k)
897+
enddo
955898

956899
! Computing probabilities
957900
if (lnRE.and.(k_abs_tot > tiny_dp)) then
@@ -970,30 +913,18 @@ subroutine opacite(lambda, p_lambda, no_scatt)
970913
enddo
971914
endif ! letape_th and not onlyLTE
972915

973-
974916
! proba absorption sur une taille donnée
917+
! TODO : I think this loop and array could be calculated only up to p_n_cells
918+
! (and then copied over every cell if needed)
975919
if (lRE_nLTE .and. (.not.low_mem_th_emission_nLTE)) then
976920
if (letape_th) then
977921
do icell=1, n_cells
978922
kabs_nLTE_CDF(grain_RE_nLTE_start-1,icell,lambda)=0.0
979-
if (lvariable_dust) then
980-
do k=grain_RE_nLTE_start, grain_RE_nLTE_end
981-
density=dust_density(k,icell) * nbre_grains(k)
982-
kabs_nLTE_CDF(k,icell,lambda) = kabs_nLTE_CDF(k-1,icell,lambda) + &
983-
C_abs(k,lambda) * density
984-
enddo !k
985-
else
986-
d1 = dust_density(1,icell)
987-
if (d1 > tiny_dp) then
988-
do k=grain_RE_nLTE_start, grain_RE_nLTE_end
989-
density=d1 * nbre_grains(k)
990-
kabs_nLTE_CDF(k,icell,lambda) = kabs_nLTE_CDF(k-1,icell,lambda) + &
991-
C_abs(k,lambda) * density
992-
enddo !k
993-
else
994-
kabs_nLTE_CDF(grain_RE_nLTE_start:grain_RE_nLTE_end,icell,lambda) = 0.0
995-
endif
996-
endif
923+
do k=grain_RE_nLTE_start, grain_RE_nLTE_end
924+
density=dust_density(p_k,icell) * nbre_grains(k)
925+
kabs_nLTE_CDF(k,icell,lambda) = kabs_nLTE_CDF(k-1,icell,lambda) + &
926+
C_abs(k,lambda) * density
927+
enddo !k
997928
if (kabs_nLTE_CDF(grain_RE_nLTE_end,icell,lambda) > tiny_real) then
998929
kabs_nLTE_CDF(:,icell,lambda) = kabs_nLTE_CDF(:,icell,lambda)/&
999930
kabs_nLTE_CDF(grain_RE_nLTE_end,icell,lambda)
@@ -1139,56 +1070,28 @@ subroutine calc_local_scattering_matrices(lambda, p_lambda)
11391070
endif
11401071
endif
11411072

1142-
if (lvariable_dust) then
1143-
do k=1,n_grains_tot
1144-
density=dust_density(k,icell) * nbre_grains(k)
1145-
if (aniso_method==1) then
1146-
! Moyennage matrice de mueller (long en cpu ) (le dernier indice est l'angle)
1147-
! tab_s11 est normalisee a Qsca --> facteur S_grain * density pour que
1148-
! tab_s11_pos soit normalisee a k_sca_tot
1149-
tab_s11_pos(:,icell,p_lambda) = tab_s11_pos(:,icell,p_lambda) + &
1150-
tab_s11(:,k,lambda) * S_grain(k) * density
1151-
if (lsepar_pola) then
1152-
tab_s12_o_s11_pos(:,icell,p_lambda) = tab_s12_o_s11_pos(:,icell,p_lambda) + &
1153-
tab_s12(:,k,lambda) * S_grain(k) * density
1154-
tab_s22_o_s11_pos(:,icell,p_lambda) = tab_s22_o_s11_pos(:,icell,p_lambda) + &
1155-
tab_s22(:,k,lambda) * S_grain(k) * density
1156-
tab_s33_o_s11_pos(:,icell,p_lambda) = tab_s33_o_s11_pos(:,icell,p_lambda) + &
1157-
tab_s33(:,k,lambda) * S_grain(k) * density
1158-
tab_s34_o_s11_pos(:,icell,p_lambda) = tab_s34_o_s11_pos(:,icell,p_lambda) + &
1159-
tab_s34(:,k,lambda) * S_grain(k) * density
1160-
tab_s44_o_s11_pos(:,icell,p_lambda) = tab_s44_o_s11_pos(:,icell,p_lambda) + &
1161-
tab_s44(:,k,lambda) * S_grain(k) * density
1162-
endif
1163-
endif !aniso_method
1164-
enddo !k
1165-
else
1166-
d1 = dust_density(1,icell)
1167-
if (d1 > tiny_dp) then
1168-
do k=1,n_grains_tot
1169-
density=d1 * nbre_grains(k)
1170-
if (aniso_method==1) then
1171-
! Moyennage matrice de mueller (long en cpu ) (le dernier indice est l'angle)
1172-
! tab_s11 est normalisee a Qsca --> facteur S_grain * density pour que
1173-
! tab_s11_pos soit normalisee a k_sca_tot
1174-
tab_s11_pos(:,icell,p_lambda) = tab_s11_pos(:,icell,p_lambda) + &
1175-
tab_s11(:,k,lambda) * S_grain(k) * density
1176-
if (lsepar_pola) then
1177-
tab_s12_o_s11_pos(:,icell,p_lambda) = tab_s12_o_s11_pos(:,icell,p_lambda) + &
1178-
tab_s12(:,k,lambda) * S_grain(k) * density
1179-
tab_s22_o_s11_pos(:,icell,p_lambda) = tab_s22_o_s11_pos(:,icell,p_lambda) + &
1180-
tab_s22(:,k,lambda) * S_grain(k) * density
1181-
tab_s33_o_s11_pos(:,icell,p_lambda) = tab_s33_o_s11_pos(:,icell,p_lambda) + &
1182-
tab_s33(:,k,lambda) * S_grain(k) * density
1183-
tab_s34_o_s11_pos(:,icell,p_lambda) = tab_s34_o_s11_pos(:,icell,p_lambda) + &
1184-
tab_s34(:,k,lambda) * S_grain(k) * density
1185-
tab_s44_o_s11_pos(:,icell,p_lambda) = tab_s44_o_s11_pos(:,icell,p_lambda) + &
1186-
tab_s44(:,k,lambda) * S_grain(k) * density
1187-
endif
1188-
endif !aniso_method
1189-
enddo !k
1190-
endif
1191-
endif
1073+
do k=1,n_grains_tot
1074+
density=dust_density(k,icell) * nbre_grains(k)
1075+
if (aniso_method==1) then
1076+
! Moyennage matrice de mueller (long en cpu ) (le dernier indice est l'angle)
1077+
! tab_s11 est normalisee a Qsca --> facteur S_grain * density pour que
1078+
! tab_s11_pos soit normalisee a k_sca_tot
1079+
tab_s11_pos(:,icell,p_lambda) = tab_s11_pos(:,icell,p_lambda) + &
1080+
tab_s11(:,k,lambda) * S_grain(k) * density
1081+
if (lsepar_pola) then
1082+
tab_s12_o_s11_pos(:,icell,p_lambda) = tab_s12_o_s11_pos(:,icell,p_lambda) + &
1083+
tab_s12(:,k,lambda) * S_grain(k) * density
1084+
tab_s22_o_s11_pos(:,icell,p_lambda) = tab_s22_o_s11_pos(:,icell,p_lambda) + &
1085+
tab_s22(:,k,lambda) * S_grain(k) * density
1086+
tab_s33_o_s11_pos(:,icell,p_lambda) = tab_s33_o_s11_pos(:,icell,p_lambda) + &
1087+
tab_s33(:,k,lambda) * S_grain(k) * density
1088+
tab_s34_o_s11_pos(:,icell,p_lambda) = tab_s34_o_s11_pos(:,icell,p_lambda) + &
1089+
tab_s34(:,k,lambda) * S_grain(k) * density
1090+
tab_s44_o_s11_pos(:,icell,p_lambda) = tab_s44_o_s11_pos(:,icell,p_lambda) + &
1091+
tab_s44(:,k,lambda) * S_grain(k) * density
1092+
endif
1093+
endif !aniso_method
1094+
enddo !k
11921095

11931096
k_sca_tot = kappa(icell,lambda) * tab_albedo_pos(icell,lambda) / fact ! We renormalize to remove the factor from opacity
11941097

@@ -1379,41 +1282,19 @@ integer function select_scattering_grain(lambda,icell, aleat) result(k)
13791282
if (aleat < 0.5) then ! We start from first grain
13801283
prob = aleat * norm
13811284
CDF = 0.0
1382-
if (lvariable_dust) then
1383-
do k=1, n_grains_tot
1384-
density = dust_density(k,icell) * nbre_grains(k)
1385-
CDF = CDF + C_sca(k,lambda) * density
1386-
if (CDF > prob) exit
1387-
enddo
1388-
else
1389-
d1 = dust_density(1,icell)
1390-
if (d1 > tiny_dp) then
1391-
do k=1, n_grains_tot
1392-
density = d1 * nbre_grains(k)
1393-
CDF = CDF + C_sca(k,lambda) * density
1394-
if (CDF > prob) exit
1395-
enddo
1396-
endif
1397-
endif
1285+
do k=1, n_grains_tot
1286+
density = dust_density(k,icell) * nbre_grains(k)
1287+
CDF = CDF + C_sca(k,lambda) * density
1288+
if (CDF > prob) exit
1289+
enddo
13981290
else ! We start from the end of the grain size distribution
13991291
prob = (1.0-aleat) * norm
14001292
CDF = 0.0
1401-
if (lvariable_dust) then
1402-
do k=n_grains_tot, 1, -1
1403-
density = dust_density(k,icell) * nbre_grains(k)
1404-
CDF = CDF + C_sca(k,lambda) * density
1405-
if (CDF > prob) exit
1406-
enddo
1407-
else
1408-
d1 = dust_density(1,icell)
1409-
if (d1 > tiny_dp) then
1410-
do k=n_grains_tot, 1, -1
1411-
density = d1 * nbre_grains(k)
1412-
CDF = CDF + C_sca(k,lambda) * density
1413-
if (CDF > prob) exit
1414-
enddo
1415-
endif
1416-
endif
1293+
do k=n_grains_tot, 1, -1
1294+
density = dust_density(k,icell) * nbre_grains(k)
1295+
CDF = CDF + C_sca(k,lambda) * density
1296+
if (CDF > prob) exit
1297+
enddo
14171298
endif
14181299
else
14191300
k = select_grainsize_high_mem(lambda,aleat, icell) ! is this ever used ?

0 commit comments

Comments
 (0)