Index: /LMDZ5/trunk/libf/dynphy_lonlat/phylmd/etat0phys_netcdf.F90
===================================================================
--- /LMDZ5/trunk/libf/dynphy_lonlat/phylmd/etat0phys_netcdf.F90	(revision 2643)
+++ /LMDZ5/trunk/libf/dynphy_lonlat/phylmd/etat0phys_netcdf.F90	(revision 2644)
@@ -111,4 +111,5 @@
   INTEGER :: flag_aerosol
   INTEGER :: flag_aerosol_strat
+  LOGICAL :: flag_bc_internal_mixture
   LOGICAL :: new_aod
   REAL    :: bl95_b0, bl95_b1
@@ -131,5 +132,5 @@
                    ok_ade, ok_aie, ok_cdnc, aerosol_couple,             &
                    flag_aerosol, flag_aerosol_strat, new_aod,           &
-                   bl95_b0, bl95_b1,                                    &
+                   flag_bc_internal_mixture, bl95_b0, bl95_b1,          &
                    read_climoz,                                         &
                    alp_offset)
Index: /LMDZ5/trunk/libf/phylmd/conf_phys_m.F90
===================================================================
--- /LMDZ5/trunk/libf/phylmd/conf_phys_m.F90	(revision 2643)
+++ /LMDZ5/trunk/libf/phylmd/conf_phys_m.F90	(revision 2644)
@@ -19,5 +19,5 @@
        ok_ade, ok_aie, ok_cdnc, aerosol_couple, &
        flag_aerosol, flag_aerosol_strat, new_aod, &
-       bl95_b0, bl95_b1,&
+       flag_bc_internal_mixture, bl95_b0, bl95_b1,&
        read_climoz, &
        alp_offset)
@@ -64,4 +64,5 @@
     ! ok_cdnc, ok cloud droplet number concentration
     ! flag_aerosol_strat : flag pour les aerosols stratos
+    ! flag_bc_internal_mixture : use BC internal mixture if true
     ! bl95_b*: parameters in the formula to link CDNC to aerosol mass conc 
     !
@@ -77,4 +78,5 @@
     INTEGER              :: flag_aerosol
     INTEGER              :: flag_aerosol_strat
+    LOGICAL              :: flag_bc_internal_mixture
     LOGICAL              :: new_aod
     REAL                 :: bl95_b0, bl95_b1
@@ -95,4 +97,5 @@
     INTEGER, SAVE       :: flag_aerosol_omp
     INTEGER, SAVE       :: flag_aerosol_strat_omp
+    LOGICAL, SAVE       :: flag_bc_internal_mixture_omp
     LOGICAL, SAVE       :: new_aod_omp
     REAL,SAVE           :: bl95_b0_omp, bl95_b1_omp
@@ -398,4 +401,15 @@
     flag_aerosol_strat_omp = 0
     CALL getin('flag_aerosol_strat',flag_aerosol_strat_omp)
+
+    !
+    !Config Key  = flag_bc_internal_mixture
+    !Config Desc = state of mixture for BC aerosols
+    ! - n = external mixture
+    ! - y = internal mixture
+    !Config Def  = n
+    !Config Help = Used in physiq.F / aeropt
+    !
+    flag_bc_internal_mixture_omp = .false.
+    CALL getin('flag_bc_internal_mixture',flag_bc_internal_mixture_omp)
 
     ! Temporary variable for testing purpose!
@@ -2136,4 +2150,5 @@
     flag_aerosol=flag_aerosol_omp
     flag_aerosol_strat=flag_aerosol_strat_omp
+    flag_bc_internal_mixture=flag_bc_internal_mixture_omp
     new_aod=new_aod_omp
     aer_type = aer_type_omp
@@ -2311,4 +2326,15 @@
     IF (ok_aie .AND. .NOT. ok_cdnc) THEN
        CALL abort_physic('conf_phys', 'ok_cdnc must be set to y if ok_aie is activated',1)
+    ENDIF
+
+    ! BC internal mixture is only possible with RRTM & NSW=6 & flag_aerosol=6 or aerosol_couple
+    IF (flag_bc_internal_mixture .AND. NSW.NE.6) THEN 
+       CALL abort_physic('conf_phys', 'flag_bc_internal_mixture can only be activated with NSW=6',1)
+    ENDIF
+    IF (flag_bc_internal_mixture .AND. iflag_rrtm.NE.1) THEN 
+       CALL abort_physic('conf_phys', 'flag_bc_internal_mixture can only be activated with RRTM',1)
+    ENDIF
+    IF (flag_bc_internal_mixture .AND. flag_aerosol.NE.6) THEN 
+       CALL abort_physic('conf_phys', 'flag_bc_internal_mixture can only be activated with flag_aerosol=6',1)
     ENDIF
 
Index: /LMDZ5/trunk/libf/phylmd/physiq_mod.F90
===================================================================
--- /LMDZ5/trunk/libf/phylmd/physiq_mod.F90	(revision 2643)
+++ /LMDZ5/trunk/libf/phylmd/physiq_mod.F90	(revision 2644)
@@ -983,4 +983,6 @@
     INTEGER, SAVE :: flag_aerosol 
     !$OMP THREADPRIVATE(flag_aerosol) 
+    LOGICAL, SAVE :: flag_bc_internal_mixture
+    !$OMP THREADPRIVATE(flag_bc_internal_mixture)
     LOGICAL, SAVE :: new_aod
     !$OMP THREADPRIVATE(new_aod) 
@@ -1140,5 +1142,5 @@
             ok_ade, ok_aie, ok_cdnc, aerosol_couple,  &
             flag_aerosol, flag_aerosol_strat, new_aod, &
-            bl95_b0, bl95_b1, &
+            flag_bc_internal_mixture, bl95_b0, bl95_b1, &
                                 ! nv flags pour la convection et les
                                 ! poches froides
@@ -3388,5 +3390,5 @@
                    !
                    CALL readaerosol_optic_rrtm( debut, aerosol_couple, &
-                        new_aod, flag_aerosol, itap, jD_cur-jD_ref, &
+                        new_aod, flag_aerosol, flag_bc_internal_mixture, itap, jD_cur-jD_ref, &
                         pdtphys, pplay, paprs, t_seri, rhcl, presnivs,  &
                         tr_seri, mass_solu_aero, mass_solu_aero_pi,  &
Index: /LMDZ5/trunk/libf/phylmd/rrtm/aeropt_5wv_rrtm.F90
===================================================================
--- /LMDZ5/trunk/libf/phylmd/rrtm/aeropt_5wv_rrtm.F90	(revision 2643)
+++ /LMDZ5/trunk/libf/phylmd/rrtm/aeropt_5wv_rrtm.F90	(revision 2644)
@@ -3,8 +3,9 @@
 !
 
-SUBROUTINE AEROPT_5WV_RRTM(&
-   pdel, m_allaer, delt, &
-   RHcl, ai, flag_aerosol, &
-   pplay, t_seri, &
+SUBROUTINE AEROPT_5WV_RRTM(  &
+   pdel, m_allaer,           &
+   RHcl, ai, flag_aerosol,   &
+   flag_bc_internal_mixture, & 
+   pplay, t_seri,            &
    tausum, tau )
 
@@ -54,8 +55,8 @@
   !
   REAL, DIMENSION(klon,klev), INTENT(in)   :: pdel
-  REAL, INTENT(in)                         :: delt
   REAL, DIMENSION(klon,klev,naero_tot), INTENT(in) :: m_allaer
   REAL, DIMENSION(klon,klev), INTENT(in)   :: RHcl     ! humidite relative ciel clair
   INTEGER,INTENT(in)                       :: flag_aerosol
+  LOGICAL,INTENT(in)                       :: flag_bc_internal_mixture
   REAL, DIMENSION(klon,klev), INTENT(in)   :: pplay
   REAL, DIMENSION(klon,klev), INTENT(in)   :: t_seri
@@ -106,8 +107,20 @@
 
   !
+  ! BC internal mixture
+  !
+  INTEGER, PARAMETER ::  nbclassbc = 5  ! Added by Rong Wang/OB for the 5 fractions
+                                       ! of BC in the soluble mode:
+                                       ! bc_content/0.001, 0.01, 0.02, 0.05, ! 0.1/
+  ! for Maxwell-Garnet internal mixture
+  ! Detailed theory can be found in R. Wang Estimation of global black carbon ! direct
+  ! radiative forcing and its uncertainty constrained by observations. J.
+  ! Geophys. Res. Atmos. Added by R. Wang and OB
+  REAL :: alpha_MG_5wv(nbre_RH,las,nbclassbc)
+
+  !
   ! Proprietes optiques
   !
-  REAL :: fact_RH(nbre_RH)
-  INTEGER :: n
+  REAL :: fact_RH(nbre_RH), BC_massfra
+  INTEGER :: n, classbc
 
 ! From here on we look at the optical parameters at 5 wavelengths:  
@@ -117,80 +130,45 @@
  DATA alpha_aers_5wv/ & 
                                 ! bc soluble 
-       7.930,7.930,7.930,7.930,7.930,7.930,     & 
-       7.930,7.930,10.893,12.618,14.550,16.613, & 
-       7.658,7.658,7.658,7.658,7.658,7.658,     & 
-       7.658,7.658,10.351,11.879,13.642,15.510, & 
-       7.195,7.195,7.195,7.195,7.195,7.195,     & 
-       7.195,7.195,9.551,10.847,12.381,13.994,  & 
-       6.736,6.736,6.736,6.736,6.736,6.736,     & 
-       6.736,6.736,8.818,9.938,11.283,12.687,   & 
-       6.277,6.277,6.277,6.277,6.277,6.277,     & 
-       6.277,6.277,8.123,9.094,10.275,11.501,   & 
+       7.930,7.930,7.930,7.930,7.930,7.930,7.930,7.930,10.893,12.618,14.550,16.613, & 
+       7.658,7.658,7.658,7.658,7.658,7.658,7.658,7.658,10.351,11.879,13.642,15.510, & 
+       7.195,7.195,7.195,7.195,7.195,7.195,7.195,7.195,9.551,10.847,12.381,13.994,  & 
+       6.736,6.736,6.736,6.736,6.736,6.736,6.736,6.736,8.818,9.938,11.283,12.687,   & 
+       6.277,6.277,6.277,6.277,6.277,6.277,6.277,6.277,8.123,9.094,10.275,11.501,   & 
                                 ! pom soluble 
-       6.676,6.676,6.676,6.676,6.710,6.934,   & 
-       7.141,7.569,8.034,8.529,9.456,10.511,  & 
-       5.109,5.109,5.109,5.109,5.189,5.535,   & 
-       5.960,6.852,8.008,9.712,12.897,19.676, & 
-       3.718,3.718,3.718,3.718,3.779,4.042,   & 
-       4.364,5.052,5.956,7.314,9.896,15.688,  & 
-       2.849,2.849,2.849,2.849,2.897,3.107,   & 
-       3.365,3.916,4.649,5.760,7.900,12.863,  & 
-       2.229,2.229,2.229,2.229,2.268,2.437,   & 
-       2.645,3.095,3.692,4.608,6.391,10.633,  & 
+       6.676,6.676,6.676,6.676,6.710,6.934,7.141,7.569,8.034,8.529,9.456,10.511,  & 
+       5.109,5.109,5.109,5.109,5.189,5.535,5.960,6.852,8.008,9.712,12.897,19.676, & 
+       3.718,3.718,3.718,3.718,3.779,4.042,4.364,5.052,5.956,7.314,9.896,15.688,  & 
+       2.849,2.849,2.849,2.849,2.897,3.107,3.365,3.916,4.649,5.760,7.900,12.863,  & 
+       2.229,2.229,2.229,2.229,2.268,2.437,2.645,3.095,3.692,4.608,6.391,10.633,  & 
                                 ! Sulfate (Accumulation) 
-       5.751,6.215,6.690,7.024,7.599,8.195,      & 
-       9.156,10.355,12.660,14.823,18.908,24.508, & 
-       4.320,4.675,5.052,5.375,5.787,6.274,      & 
-       7.066,8.083,10.088,12.003,15.697,21.133,  & 
-       3.079,3.351,3.639,3.886,4.205,4.584,      & 
-       5.206,6.019,7.648,9.234,12.391,17.220,    & 
-       2.336,2.552,2.781,2.979,3.236,3.540,      & 
-       4.046,4.711,6.056,7.388,10.093,14.313,    & 
-       1.777,1.949,2.134,2.292,2.503,2.751,      & 
-       3.166,3.712,4.828,5.949,8.264,11.922,     & 
+       5.751,6.215,6.690,7.024,7.599,8.195,9.156,10.355,12.660,14.823,18.908,24.508, & 
+       4.320,4.675,5.052,5.375,5.787,6.274,7.066,8.083,10.088,12.003,15.697,21.133,  & 
+       3.079,3.351,3.639,3.886,4.205,4.584,5.206,6.019,7.648,9.234,12.391,17.220,    & 
+       2.336,2.552,2.781,2.979,3.236,3.540,4.046,4.711,6.056,7.388,10.093,14.313,    & 
+       1.777,1.949,2.134,2.292,2.503,2.751,3.166,3.712,4.828,5.949,8.264,11.922,     & 
                                 ! Sulfate (Coarse) 
-       5.751,6.215,6.690,7.024,7.599,8.195,      & 
-       9.156,10.355,12.660,14.823,18.908,24.508, & 
-       4.320,4.675,5.052,5.375,5.787,6.274,      & 
-       7.066,8.083,10.088,12.003,15.697,21.133,  & 
-       3.079,3.351,3.639,3.886,4.205,4.584,      & 
-       5.206,6.019,7.648,9.234,12.391,17.220,    & 
-       2.336,2.552,2.781,2.979,3.236,3.540,      & 
-       4.046,4.711,6.056,7.388,10.093,14.313,    & 
-       1.777,1.949,2.134,2.292,2.503,2.751,      & 
-       3.166,3.712,4.828,5.949,8.264,11.922,     & 
+       5.751,6.215,6.690,7.024,7.599,8.195,9.156,10.355,12.660,14.823,18.908,24.508, & 
+       4.320,4.675,5.052,5.375,5.787,6.274,7.066,8.083,10.088,12.003,15.697,21.133,  & 
+       3.079,3.351,3.639,3.886,4.205,4.584,5.206,6.019,7.648,9.234,12.391,17.220,    & 
+       2.336,2.552,2.781,2.979,3.236,3.540,4.046,4.711,6.056,7.388,10.093,14.313,    & 
+       1.777,1.949,2.134,2.292,2.503,2.751,3.166,3.712,4.828,5.949,8.264,11.922,     & 
                            ! seasalt seasalt Super Coarse Soluble (SS)
-        0.218, 0.272, 0.293, 0.316, 0.343, 0.380, &
-        0.429, 0.501, 0.636, 0.755, 0.967, 1.495, & 
-        0.221, 0.275, 0.297, 0.320, 0.348, 0.383, &
-        0.432, 0.509, 0.640, 0.759, 0.972, 1.510, &
-        0.224, 0.279, 0.301, 0.324, 0.352, 0.388, &
-        0.438, 0.514, 0.647, 0.768, 0.985, 1.514, &
-        0.227, 0.282, 0.303, 0.327, 0.356, 0.392, &
-        0.441, 0.518, 0.652, 0.770, 0.987, 1.529, &
-        0.230, 0.285, 0.306, 0.330, 0.359, 0.396, &
-        0.446, 0.522, 0.656, 0.777, 0.993, 1.539, &
+        0.218, 0.272, 0.293, 0.316, 0.343, 0.380, 0.429, 0.501, 0.636, 0.755, 0.967, 1.495, & 
+        0.221, 0.275, 0.297, 0.320, 0.348, 0.383, 0.432, 0.509, 0.640, 0.759, 0.972, 1.510, &
+        0.224, 0.279, 0.301, 0.324, 0.352, 0.388, 0.438, 0.514, 0.647, 0.768, 0.985, 1.514, &
+        0.227, 0.282, 0.303, 0.327, 0.356, 0.392, 0.441, 0.518, 0.652, 0.770, 0.987, 1.529, &
+        0.230, 0.285, 0.306, 0.330, 0.359, 0.396, 0.446, 0.522, 0.656, 0.777, 0.993, 1.539, &
                            ! seasalt seasalt Coarse Soluble (CS)      
-        0.578, 0.706, 0.756, 0.809, 0.876, 0.964, &
-        1.081, 1.256, 1.577, 1.858, 2.366, 3.613, &
-        0.598, 0.725, 0.779, 0.833, 0.898, 0.990, &
-        1.109, 1.290, 1.609, 1.889, 2.398, 3.682, &
-        0.619, 0.750, 0.802, 0.857, 0.927, 1.022, &
-        1.141, 1.328, 1.648, 1.939, 2.455, 3.729, &
-        0.633, 0.767, 0.820, 0.879, 0.948, 1.044, &
-        1.167, 1.353, 1.683, 1.969, 2.491, 3.785, & 
-        0.648, 0.785, 0.838, 0.896, 0.967, 1.066, &
-        1.192, 1.381, 1.714, 2.006, 2.531, 3.836, &
+        0.578, 0.706, 0.756, 0.809, 0.876, 0.964, 1.081, 1.256, 1.577, 1.858, 2.366, 3.613, &
+        0.598, 0.725, 0.779, 0.833, 0.898, 0.990, 1.109, 1.290, 1.609, 1.889, 2.398, 3.682, &
+        0.619, 0.750, 0.802, 0.857, 0.927, 1.022, 1.141, 1.328, 1.648, 1.939, 2.455, 3.729, &
+        0.633, 0.767, 0.820, 0.879, 0.948, 1.044, 1.167, 1.353, 1.683, 1.969, 2.491, 3.785, & 
+        0.648, 0.785, 0.838, 0.896, 0.967, 1.066, 1.192, 1.381, 1.714, 2.006, 2.531, 3.836, &
                            ! seasalt seasalt Accumulation Soluble (AS)
-        4.432, 5.899, 6.505, 7.166, 7.964, 7.962, &
-        9.232,11.257,14.979,18.337,24.223,37.811, &
-        3.298, 4.569, 5.110, 5.709, 6.446, 6.268, &
-        7.396, 9.246,12.787,16.113,22.197,37.136, &
-        2.340, 3.358, 3.803, 4.303, 4.928, 4.696, &
-        5.629, 7.198,10.308,13.342,19.120,34.296, &
-        1.789, 2.626, 2.999, 3.422, 3.955, 3.730, &
-        4.519, 5.864, 8.593,11.319,16.653,31.331, &
-        1.359, 2.037, 2.343, 2.693, 3.139, 2.940, &
-        3.596, 4.729, 7.076, 9.469,14.266,28.043 /
+        4.432, 5.899, 6.505, 7.166, 7.964, 7.962, 9.232,11.257,14.979,18.337,24.223,37.811, &
+        3.298, 4.569, 5.110, 5.709, 6.446, 6.268, 7.396, 9.246,12.787,16.113,22.197,37.136, &
+        2.340, 3.358, 3.803, 4.303, 4.928, 4.696, 5.629, 7.198,10.308,13.342,19.120,34.296, &
+        1.789, 2.626, 2.999, 3.422, 3.955, 3.730, 4.519, 5.864, 8.593,11.319,16.653,31.331, &
+        1.359, 2.037, 2.343, 2.693, 3.139, 2.940, 3.596, 4.729, 7.076, 9.469,14.266,28.043 /
 
   DATA alpha_aeri_5wv/ &
@@ -201,4 +179,38 @@
                                  ! pom insoluble 
         5.042, 3.101, 1.890, 1.294, 0.934/ 
+
+  ! internal mixture
+  DATA alpha_MG_5wv/ &
+     ! bc content = 0.001
+       7.930,7.930,7.930,7.930,7.930,7.930,7.930,7.930,10.893,12.618,14.550,16.613, & 
+       7.658,7.658,7.658,7.658,7.658,7.658,7.658,7.658,10.351,11.879,13.642,15.510, & 
+       7.195,7.195,7.195,7.195,7.195,7.195,7.195,7.195,9.551,10.847,12.381,13.994,  & 
+       6.736,6.736,6.736,6.736,6.736,6.736,6.736,6.736,8.818,9.938,11.283,12.687,   & 
+       6.277,6.277,6.277,6.277,6.277,6.277,6.277,6.277,8.123,9.094,10.275,11.501,   & 
+     ! bc content = 0.01
+       7.930,7.930,7.930,7.930,7.930,7.930,7.930,7.930,10.893,12.618,14.550,16.613, & 
+       7.658,7.658,7.658,7.658,7.658,7.658,7.658,7.658,10.351,11.879,13.642,15.510, & 
+       7.195,7.195,7.195,7.195,7.195,7.195,7.195,7.195,9.551,10.847,12.381,13.994,  & 
+       6.736,6.736,6.736,6.736,6.736,6.736,6.736,6.736,8.818,9.938,11.283,12.687,   & 
+       6.277,6.277,6.277,6.277,6.277,6.277,6.277,6.277,8.123,9.094,10.275,11.501,   & 
+     ! bc content = 0.02
+       7.930,7.930,7.930,7.930,7.930,7.930,7.930,7.930,10.893,12.618,14.550,16.613, & 
+       7.658,7.658,7.658,7.658,7.658,7.658,7.658,7.658,10.351,11.879,13.642,15.510, & 
+       7.195,7.195,7.195,7.195,7.195,7.195,7.195,7.195,9.551,10.847,12.381,13.994,  & 
+       6.736,6.736,6.736,6.736,6.736,6.736,6.736,6.736,8.818,9.938,11.283,12.687,   & 
+       6.277,6.277,6.277,6.277,6.277,6.277,6.277,6.277,8.123,9.094,10.275,11.501,   & 
+     ! bc content = 0.05
+       7.930,7.930,7.930,7.930,7.930,7.930,7.930,7.930,10.893,12.618,14.550,16.613, & 
+       7.658,7.658,7.658,7.658,7.658,7.658,7.658,7.658,10.351,11.879,13.642,15.510, & 
+       7.195,7.195,7.195,7.195,7.195,7.195,7.195,7.195,9.551,10.847,12.381,13.994,  & 
+       6.736,6.736,6.736,6.736,6.736,6.736,6.736,6.736,8.818,9.938,11.283,12.687,   & 
+       6.277,6.277,6.277,6.277,6.277,6.277,6.277,6.277,8.123,9.094,10.275,11.501,   & 
+     ! bc content = 0.10
+       7.930,7.930,7.930,7.930,7.930,7.930,7.930,7.930,10.893,12.618,14.550,16.613, & 
+       7.658,7.658,7.658,7.658,7.658,7.658,7.658,7.658,10.351,11.879,13.642,15.510, & 
+       7.195,7.195,7.195,7.195,7.195,7.195,7.195,7.195,9.551,10.847,12.381,13.994,  & 
+       6.736,6.736,6.736,6.736,6.736,6.736,6.736,6.736,8.818,9.938,11.283,12.687,   & 
+       6.277,6.277,6.277,6.277,6.277,6.277,6.277,6.277,8.123,9.094,10.275,11.501 /
+
   ! 
   ! Initialisations
@@ -251,5 +263,5 @@
      aerosol_name(8) = id_ASSSM_phy
      aerosol_name(9) = id_CIDUSTM_phy
-     aerosol_name(10) = id_CSSO4M_phy
+     aerosol_name(10)= id_CSSO4M_phy
   ENDIF
 
@@ -289,9 +301,9 @@
         soluble=.TRUE.
         spsol=3
-        fac=1.375    ! (NH4)2-SO4/SO4 132/96 mass conversion factor for OD
+        fac=1.375    ! (NH4)2-SO4/SO4 132/96 mass conversion factor for AOD
     ELSEIF (aerosol_name(m).EQ.id_CSSO4M_phy) THEN
         soluble=.TRUE.
         spsol=4
-        fac=1.375    ! (NH4)2-SO4/SO4 132/96 mass conversion factor for OD
+        fac=1.375    ! (NH4)2-SO4/SO4 132/96 mass conversion factor for AOD
     ELSEIF (aerosol_name(m).EQ.id_SSSSM_phy) THEN 
         soluble=.TRUE.
@@ -325,4 +337,35 @@
       IF (soluble) THEN            ! For soluble aerosol
 
+        !--treat special case of soluble BC internal mixture
+        IF (spsol.EQ.1 .AND. flag_bc_internal_mixture) THEN
+
+          DO k=1, klev
+            DO i=1, klon
+
+             BC_massfra = m_allaer(i,k,id_ASBCM_phy)/(m_allaer(i,k,id_ASBCM_phy)+m_allaer(i,k,id_ASSO4M_phy))
+
+             IF (BC_massfra.GE.0.10) THEN
+               classbc = 5
+             ELSEIF  (BC_massfra.GE.0.05) THEN
+               classbc = 4
+             ELSEIF  (BC_massfra.GE.0.02) THEN
+               classbc = 3
+             ELSEIF  (BC_massfra.GE.0.01) THEN
+               classbc = 2
+             ELSE
+               classbc = 1
+             ENDIF
+
+              tau_ae5wv_int = alpha_MG_5wv(RH_num(i,k),la,classbc)+DELTA(i,k)* &
+                             (alpha_MG_5wv(RH_num(i,k)+1,la,classbc) - & 
+                              alpha_MG_5wv(RH_num(i,k),la,classbc))
+              tau(i,k,la,aerindex) = m_allaer(i,k,aerindex)/1.e6*zdh(i,k)*tau_ae5wv_int*fac
+              tausum(i,la,aerindex)=tausum(i,la,aerindex)+tau(i,k,la,aerindex)
+            ENDDO
+          ENDDO
+
+        !--other cases of soluble aerosols
+        ELSE
+
           DO k=1, klev
             DO i=1, klon
@@ -334,6 +377,9 @@
             ENDDO
           ENDDO
+
+        ENDIF
  
-      ELSE                         ! For insoluble aerosol
+      ! cases of insoluble aerosol
+      ELSE                         
 
         DO k=1, klev
Index: /LMDZ5/trunk/libf/phylmd/rrtm/aeropt_6bands_rrtm.F90
===================================================================
--- /LMDZ5/trunk/libf/phylmd/rrtm/aeropt_6bands_rrtm.F90	(revision 2643)
+++ /LMDZ5/trunk/libf/phylmd/rrtm/aeropt_6bands_rrtm.F90	(revision 2644)
@@ -3,8 +3,8 @@
 !
 SUBROUTINE AEROPT_6BANDS_RRTM ( &
-     pdel, m_allaer, delt, RHcl, &
+     pdel, m_allaer, RHcl, &
      tau_allaer, piz_allaer, &
      cg_allaer, m_allaer_pi, &
-     flag_aerosol, zrho )
+     flag_aerosol, flag_bc_internal_mixture, zrho )
 
   USE dimphy
@@ -28,9 +28,9 @@
   !
   REAL, DIMENSION(klon,klev),     INTENT(in)  :: pdel
-  REAL,                           INTENT(in)  :: delt
   REAL, DIMENSION(klon,klev,naero_tot),   INTENT(in)  :: m_allaer
   REAL, DIMENSION(klon,klev,naero_tot),   INTENT(in)  :: m_allaer_pi
   REAL, DIMENSION(klon,klev),     INTENT(in)  :: RHcl       ! humidite relative ciel clair
   INTEGER,                        INTENT(in)  :: flag_aerosol
+  LOGICAL,                        INTENT(in)  :: flag_bc_internal_mixture
   REAL, DIMENSION(klon,klev),     INTENT(in)  :: zrho
   !
@@ -71,7 +71,10 @@
 
   REAL, DIMENSION(klon,klev,naero_tot,nbands_sw_rrtm) ::  tau_ae
-  REAL, DIMENSION(klon,klev,naero_tot,nbands_sw_rrtm) ::  tau_ae_pi
   REAL, DIMENSION(klon,klev,naero_tot,nbands_sw_rrtm) ::  piz_ae
   REAL, DIMENSION(klon,klev,naero_tot,nbands_sw_rrtm) ::  cg_ae
+
+  REAL, DIMENSION(klon,klev,naero_tot,nbands_sw_rrtm) ::  tau_ae_pi
+  REAL, DIMENSION(klon,klev,id_ASBCM_phy:id_ASBCM_phy,nbands_sw_rrtm) :: piz_ae_pi
+  REAL, DIMENSION(klon,klev,id_ASBCM_phy:id_ASBCM_phy,nbands_sw_rrtm) :: cg_ae_pi
   !
   ! Proprietes optiques
@@ -83,7 +86,24 @@
   REAL:: piz_aers_6bands(nbre_RH,nbands_sw_rrtm,naero_soluble)     !-- unit
   REAL:: piz_aeri_6bands(nbands_sw_rrtm,naero_insoluble)        !-- unit
-
-  INTEGER :: id
-  REAL :: tmp_var, tmp_var_pi
+  ! 
+  ! BC internal mixture
+  !
+  INTEGER, PARAMETER ::  nbclassbc = 5  ! Added by Rong Wang/OB for the 5 fractions
+                                       ! of BC in the soluble mode:
+                                       ! bc_content/0.001, 0.01, 0.02, 0.05, 0.1/
+  ! for Maxwell-Garnet internal mixture
+  ! Detailed theory can be found in R. Wang Estimation of global black carbon direct
+  ! radiative forcing and its uncertainty constrained by observations. J.
+  ! Geophys. Res. Atmos. Added by R. Wang and OB
+  REAL :: alpha_MG_6bands(nbre_RH,nbands_sw_rrtm,nbclassbc)
+  REAL :: cg_MG_6bands(nbre_RH,nbands_sw_rrtm,nbclassbc)
+  REAL :: piz_MG_6bands(nbre_RH,nbands_sw_rrtm,nbclassbc)
+  !
+  INTEGER :: id, classbc, classbc_pi
+  REAL :: tmp_var, tmp_var_pi, BC_massfra, BC_massfra_pi
+
+  !
+  REAL, PARAMETER :: tau_min = 1.e-15
+!  REAL, PARAMETER :: tau_min = 1.e-7
 
 !***************************************************************************
@@ -268,4 +288,118 @@
   0.973, 0.973, 0.972, 0.940, 0.816, 0.663 /
 
+! Added by R. Wang (July 31 2016)
+! properties for BC assuming Maxwell-Garnett rule and internal mixture
+
+  DATA alpha_MG_6bands/ &
+     ! bc content = 0.001
+  6.497, 6.497, 6.497, 6.497, 6.497, 7.160, 7.875, 9.356,10.811,10.974,11.149,12.734, &
+  6.497, 6.497, 6.497, 6.497, 6.497, 7.160, 7.875, 9.356,10.811,10.974,11.149,12.734, &
+  5.900, 5.900, 5.900, 5.900, 5.900, 6.502, 7.151, 8.496, 9.818, 9.965,10.124,11.564, &
+  4.284, 4.284, 4.284, 4.284, 4.284, 4.721, 5.193, 6.169, 7.129, 7.236, 7.352, 8.397, &
+  2.163, 2.163, 2.163, 2.163, 2.163, 2.384, 2.622, 3.115, 3.600, 3.654, 3.712, 4.240, &
+  0.966, 0.966, 0.966, 0.966, 0.966, 1.065, 1.171, 1.392, 1.608, 1.632, 1.658, 1.894, &
+     ! bc content = 0.01
+  6.497, 6.497, 6.497, 6.497, 6.497, 7.160, 7.875, 9.356,10.811,10.974,11.149,12.734, &
+  6.497, 6.497, 6.497, 6.497, 6.497, 7.160, 7.875, 9.356,10.811,10.974,11.149,12.734, &
+  5.900, 5.900, 5.900, 5.900, 5.900, 6.502, 7.151, 8.496, 9.818, 9.965,10.124,11.564, &
+  4.284, 4.284, 4.284, 4.284, 4.284, 4.721, 5.193, 6.169, 7.129, 7.236, 7.352, 8.397, &
+  2.163, 2.163, 2.163, 2.163, 2.163, 2.384, 2.622, 3.115, 3.600, 3.654, 3.712, 4.240, &
+  0.966, 0.966, 0.966, 0.966, 0.966, 1.065, 1.171, 1.392, 1.608, 1.632, 1.658, 1.894, &
+     ! bc content = 0.02
+  6.497, 6.497, 6.497, 6.497, 6.497, 7.160, 7.875, 9.356,10.811,10.974,11.149,12.734, &
+  6.497, 6.497, 6.497, 6.497, 6.497, 7.160, 7.875, 9.356,10.811,10.974,11.149,12.734, &
+  5.900, 5.900, 5.900, 5.900, 5.900, 6.502, 7.151, 8.496, 9.818, 9.965,10.124,11.564, &
+  4.284, 4.284, 4.284, 4.284, 4.284, 4.721, 5.193, 6.169, 7.129, 7.236, 7.352, 8.397, &
+  2.163, 2.163, 2.163, 2.163, 2.163, 2.384, 2.622, 3.115, 3.600, 3.654, 3.712, 4.240, &
+  0.966, 0.966, 0.966, 0.966, 0.966, 1.065, 1.171, 1.392, 1.608, 1.632, 1.658, 1.894, &
+     ! bc content = 0.05
+  6.497, 6.497, 6.497, 6.497, 6.497, 7.160, 7.875, 9.356,10.811,10.974,11.149,12.734, &
+  6.497, 6.497, 6.497, 6.497, 6.497, 7.160, 7.875, 9.356,10.811,10.974,11.149,12.734, &
+  5.900, 5.900, 5.900, 5.900, 5.900, 6.502, 7.151, 8.496, 9.818, 9.965,10.124,11.564, &
+  4.284, 4.284, 4.284, 4.284, 4.284, 4.721, 5.193, 6.169, 7.129, 7.236, 7.352, 8.397, &
+  2.163, 2.163, 2.163, 2.163, 2.163, 2.384, 2.622, 3.115, 3.600, 3.654, 3.712, 4.240, &
+  0.966, 0.966, 0.966, 0.966, 0.966, 1.065, 1.171, 1.392, 1.608, 1.632, 1.658, 1.894, &
+     ! bc content = 0.10
+  6.497, 6.497, 6.497, 6.497, 6.497, 7.160, 7.875, 9.356,10.811,10.974,11.149,12.734, &
+  6.497, 6.497, 6.497, 6.497, 6.497, 7.160, 7.875, 9.356,10.811,10.974,11.149,12.734, &
+  5.900, 5.900, 5.900, 5.900, 5.900, 6.502, 7.151, 8.496, 9.818, 9.965,10.124,11.564, &
+  4.284, 4.284, 4.284, 4.284, 4.284, 4.721, 5.193, 6.169, 7.129, 7.236, 7.352, 8.397, &
+  2.163, 2.163, 2.163, 2.163, 2.163, 2.384, 2.622, 3.115, 3.600, 3.654, 3.712, 4.240, &
+  0.966, 0.966, 0.966, 0.966, 0.966, 1.065, 1.171, 1.392, 1.608, 1.632, 1.658, 1.894 /
+
+  DATA cg_MG_6bands/ &
+     ! bc content = 0.001
+  0.721, 0.721, 0.721, 0.729, 0.735, 0.741, 0.746, 0.754, 0.762, 0.766, 0.769, 0.775, &
+  0.721, 0.721, 0.721, 0.729, 0.735, 0.741, 0.746, 0.754, 0.762, 0.766, 0.769, 0.775, &
+  0.643, 0.643, 0.643, 0.654, 0.662, 0.670, 0.677, 0.688, 0.698, 0.704, 0.707, 0.715, &
+  0.513, 0.513, 0.513, 0.522, 0.530, 0.536, 0.542, 0.552, 0.560, 0.565, 0.568, 0.575, &
+  0.321, 0.321, 0.321, 0.323, 0.325, 0.327, 0.328, 0.331, 0.333, 0.334, 0.335, 0.337, &
+  0.153, 0.153, 0.153, 0.149, 0.145, 0.142, 0.139, 0.135, 0.130, 0.128, 0.127, 0.123, &
+     ! bc content = 0.01
+  0.721, 0.721, 0.721, 0.729, 0.735, 0.741, 0.746, 0.754, 0.762, 0.766, 0.769, 0.775, &
+  0.721, 0.721, 0.721, 0.729, 0.735, 0.741, 0.746, 0.754, 0.762, 0.766, 0.769, 0.775, &
+  0.643, 0.643, 0.643, 0.654, 0.662, 0.670, 0.677, 0.688, 0.698, 0.704, 0.707, 0.715, &
+  0.513, 0.513, 0.513, 0.522, 0.530, 0.536, 0.542, 0.552, 0.560, 0.565, 0.568, 0.575, &
+  0.321, 0.321, 0.321, 0.323, 0.325, 0.327, 0.328, 0.331, 0.333, 0.334, 0.335, 0.337, &
+  0.153, 0.153, 0.153, 0.149, 0.145, 0.142, 0.139, 0.135, 0.130, 0.128, 0.127, 0.123, &
+     ! bc content = 0.02
+  0.721, 0.721, 0.721, 0.729, 0.735, 0.741, 0.746, 0.754, 0.762, 0.766, 0.769, 0.775, &
+  0.721, 0.721, 0.721, 0.729, 0.735, 0.741, 0.746, 0.754, 0.762, 0.766, 0.769, 0.775, &
+  0.643, 0.643, 0.643, 0.654, 0.662, 0.670, 0.677, 0.688, 0.698, 0.704, 0.707, 0.715, &
+  0.513, 0.513, 0.513, 0.522, 0.530, 0.536, 0.542, 0.552, 0.560, 0.565, 0.568, 0.575, &
+  0.321, 0.321, 0.321, 0.323, 0.325, 0.327, 0.328, 0.331, 0.333, 0.334, 0.335, 0.337, &
+  0.153, 0.153, 0.153, 0.149, 0.145, 0.142, 0.139, 0.135, 0.130, 0.128, 0.127, 0.123, &
+     ! bc content = 0.05
+  0.721, 0.721, 0.721, 0.729, 0.735, 0.741, 0.746, 0.754, 0.762, 0.766, 0.769, 0.775, &
+  0.721, 0.721, 0.721, 0.729, 0.735, 0.741, 0.746, 0.754, 0.762, 0.766, 0.769, 0.775, &
+  0.643, 0.643, 0.643, 0.654, 0.662, 0.670, 0.677, 0.688, 0.698, 0.704, 0.707, 0.715, &
+  0.513, 0.513, 0.513, 0.522, 0.530, 0.536, 0.542, 0.552, 0.560, 0.565, 0.568, 0.575, &
+  0.321, 0.321, 0.321, 0.323, 0.325, 0.327, 0.328, 0.331, 0.333, 0.334, 0.335, 0.337, &
+  0.153, 0.153, 0.153, 0.149, 0.145, 0.142, 0.139, 0.135, 0.130, 0.128, 0.127, 0.123, &
+     ! bc content = 0.10
+  0.721, 0.721, 0.721, 0.729, 0.735, 0.741, 0.746, 0.754, 0.762, 0.766, 0.769, 0.775, &
+  0.721, 0.721, 0.721, 0.729, 0.735, 0.741, 0.746, 0.754, 0.762, 0.766, 0.769, 0.775, &
+  0.643, 0.643, 0.643, 0.654, 0.662, 0.670, 0.677, 0.688, 0.698, 0.704, 0.707, 0.715, &
+  0.513, 0.513, 0.513, 0.522, 0.530, 0.536, 0.542, 0.552, 0.560, 0.565, 0.568, 0.575, &
+  0.321, 0.321, 0.321, 0.323, 0.325, 0.327, 0.328, 0.331, 0.333, 0.334, 0.335, 0.337, &
+  0.153, 0.153, 0.153, 0.149, 0.145, 0.142, 0.139, 0.135, 0.130, 0.128, 0.127, 0.123 /
+
+  DATA piz_MG_6bands/ &
+     ! bc content = 0.001
+  0.460, 0.460, 0.460, 0.460, 0.460, 0.534, 0.594, 0.688, 0.748, 0.754, 0.760, 0.803, &
+  0.460, 0.460, 0.460, 0.460, 0.460, 0.534, 0.594, 0.688, 0.748, 0.754, 0.760, 0.803, &
+  0.445, 0.445, 0.445, 0.445, 0.445, 0.521, 0.583, 0.679, 0.741, 0.747, 0.753, 0.798, &
+  0.394, 0.394, 0.394, 0.394, 0.394, 0.477, 0.545, 0.649, 0.718, 0.724, 0.730, 0.779, &
+  0.267, 0.267, 0.267, 0.267, 0.267, 0.365, 0.446, 0.571, 0.652, 0.660, 0.667, 0.725, &
+  0.121, 0.121, 0.121, 0.121, 0.121, 0.139, 0.155, 0.178, 0.193, 0.195, 0.196, 0.207, &
+     ! bc content = 0.01
+  0.460, 0.460, 0.460, 0.460, 0.460, 0.534, 0.594, 0.688, 0.748, 0.754, 0.760, 0.803, &
+  0.460, 0.460, 0.460, 0.460, 0.460, 0.534, 0.594, 0.688, 0.748, 0.754, 0.760, 0.803, &
+  0.445, 0.445, 0.445, 0.445, 0.445, 0.521, 0.583, 0.679, 0.741, 0.747, 0.753, 0.798, &
+  0.394, 0.394, 0.394, 0.394, 0.394, 0.477, 0.545, 0.649, 0.718, 0.724, 0.730, 0.779, &
+  0.267, 0.267, 0.267, 0.267, 0.267, 0.365, 0.446, 0.571, 0.652, 0.660, 0.667, 0.725, &
+  0.121, 0.121, 0.121, 0.121, 0.121, 0.139, 0.155, 0.178, 0.193, 0.195, 0.196, 0.207, &
+     ! bc content = 0.02
+  0.460, 0.460, 0.460, 0.460, 0.460, 0.534, 0.594, 0.688, 0.748, 0.754, 0.760, 0.803, &
+  0.460, 0.460, 0.460, 0.460, 0.460, 0.534, 0.594, 0.688, 0.748, 0.754, 0.760, 0.803, &
+  0.445, 0.445, 0.445, 0.445, 0.445, 0.521, 0.583, 0.679, 0.741, 0.747, 0.753, 0.798, &
+  0.394, 0.394, 0.394, 0.394, 0.394, 0.477, 0.545, 0.649, 0.718, 0.724, 0.730, 0.779, &
+  0.267, 0.267, 0.267, 0.267, 0.267, 0.365, 0.446, 0.571, 0.652, 0.660, 0.667, 0.725, &
+  0.121, 0.121, 0.121, 0.121, 0.121, 0.139, 0.155, 0.178, 0.193, 0.195, 0.196, 0.207, &
+     ! bc content = 0.05
+  0.460, 0.460, 0.460, 0.460, 0.460, 0.534, 0.594, 0.688, 0.748, 0.754, 0.760, 0.803, &
+  0.460, 0.460, 0.460, 0.460, 0.460, 0.534, 0.594, 0.688, 0.748, 0.754, 0.760, 0.803, &
+  0.445, 0.445, 0.445, 0.445, 0.445, 0.521, 0.583, 0.679, 0.741, 0.747, 0.753, 0.798, &
+  0.394, 0.394, 0.394, 0.394, 0.394, 0.477, 0.545, 0.649, 0.718, 0.724, 0.730, 0.779, &
+  0.267, 0.267, 0.267, 0.267, 0.267, 0.365, 0.446, 0.571, 0.652, 0.660, 0.667, 0.725, &
+  0.121, 0.121, 0.121, 0.121, 0.121, 0.139, 0.155, 0.178, 0.193, 0.195, 0.196, 0.207, &
+     ! bc content = 0.10
+  0.460, 0.460, 0.460, 0.460, 0.460, 0.534, 0.594, 0.688, 0.748, 0.754, 0.760, 0.803, &
+  0.460, 0.460, 0.460, 0.460, 0.460, 0.534, 0.594, 0.688, 0.748, 0.754, 0.760, 0.803, &
+  0.445, 0.445, 0.445, 0.445, 0.445, 0.521, 0.583, 0.679, 0.741, 0.747, 0.753, 0.798, &
+  0.394, 0.394, 0.394, 0.394, 0.394, 0.477, 0.545, 0.649, 0.718, 0.724, 0.730, 0.779, &
+  0.267, 0.267, 0.267, 0.267, 0.267, 0.365, 0.446, 0.571, 0.652, 0.660, 0.667, 0.725, &
+  0.121, 0.121, 0.121, 0.121, 0.121, 0.139, 0.155, 0.178, 0.193, 0.195, 0.196, 0.207 /
+
 !----BEGINNING OF CALCULATIONS
 
@@ -342,7 +476,10 @@
 
   tau_ae(:,:,:,:)=0.
-  tau_ae_pi(:,:,:,:)=0.
   piz_ae(:,:,:,:)=0.
   cg_ae(:,:,:,:)=0.
+
+  tau_ae_pi(:,:,:,:)=0.
+  piz_ae_pi(:,:,:,:)=0.
+  cg_ae_pi(:,:,:,:)=0.
     
   DO m=1,nb_aer   ! tau is only computed for each mass
@@ -357,9 +494,9 @@
         soluble=.TRUE.
         spsol=3
-        fac=1.375    ! (NH4)2-SO4/SO4 132/96 mass conversion factor for OD
+        fac=1.375    ! (NH4)2-SO4/SO4 132/96 mass conversion factor for AOD
      ELSEIF  (aerosol_name(m).EQ.id_CSSO4M_phy) THEN
         soluble=.TRUE.
         spsol=4
-        fac=1.375    ! (NH4)2-SO4/SO4 132/96 mass conversion factor for OD
+        fac=1.375    ! (NH4)2-SO4/SO4 132/96 mass conversion factor for AOD
      ELSEIF (aerosol_name(m).EQ.id_SSSSM_phy) THEN 
          soluble=.TRUE.
@@ -384,35 +521,126 @@
      ENDIF
 
+    !--shortname for aerosol index
     id=aerosol_name(m)
 
     IF (soluble) THEN
 
-       DO k=1, klev
-         DO i=1, klon
-           tmp_var=m_allaer(i,k,spsol)/1.e6*zdh(i,k)*fac
-           tmp_var_pi=m_allaer_pi(i,k,spsol)/1.e6*zdh(i,k)*fac
-
-           DO inu=1,NSW
-
-             tau_ae2b_int= alpha_aers_6bands(RH_num(i,k),inu,spsol)+ & 
-                           delta(i,k)* (alpha_aers_6bands(RH_num(i,k)+1,inu,spsol) - & 
-                           alpha_aers_6bands(RH_num(i,k),inu,spsol))
+       !--here we treat the special case of soluble BC internal mixture with Maxwell-Garnett rule
+       IF (spsol.EQ.1 .AND. flag_bc_internal_mixture) THEN
+
+         DO k=1, klev
+           DO i=1, klon
+
+             tmp_var=m_allaer(i,k,spsol)/1.e6*zdh(i,k)*fac
+             tmp_var_pi=m_allaer_pi(i,k,spsol)/1.e6*zdh(i,k)*fac
+
+             ! Calculate the dry BC/(BC+SUL) mass ratio for all (natural+anthropogenic) aerosols
+             BC_massfra = m_allaer(i,k,id_ASBCM_phy)/(m_allaer(i,k,id_ASBCM_phy)+m_allaer(i,k,id_ASSO4M_phy))
+
+             IF (BC_massfra.GE.0.10) THEN
+               classbc = 5
+             ELSEIF  (BC_massfra.GE.0.05) THEN
+               classbc = 4
+             ELSEIF  (BC_massfra.GE.0.02) THEN
+               classbc = 3
+             ELSEIF  (BC_massfra.GE.0.01) THEN
+               classbc = 2
+             ELSE
+               classbc = 1
+             ENDIF
+
+             ! Calculate the dry BC/(BC+SUL) mass ratio for natural aerosols
+             BC_massfra_pi = m_allaer_pi(i,k,id_ASBCM_phy)/(m_allaer_pi(i,k,id_ASBCM_phy)+m_allaer_pi(i,k,id_ASSO4M_phy))
+
+             IF (BC_massfra_pi.GE.0.10) THEN
+               classbc_pi = 5
+             ELSEIF  (BC_massfra_pi.GE.0.05) THEN
+               classbc_pi = 4
+             ELSEIF  (BC_massfra_pi.GE.0.02) THEN
+               classbc_pi = 3
+             ELSEIF  (BC_massfra_pi.GE.0.01) THEN
+               classbc_pi = 2
+             ELSE
+               classbc_pi = 1
+             ENDIF
+
+             DO inu=1,NSW
+
+               !--all aerosols
+               tau_ae2b_int= alpha_MG_6bands(RH_num(i,k),inu,classbc)+                 &
+                             delta(i,k)* (alpha_MG_6bands(RH_num(i,k)+1,inu,classbc) - & 
+                             alpha_MG_6bands(RH_num(i,k),inu,classbc))
                    
-             piz_ae2b_int = piz_aers_6bands(RH_num(i,k),inu,spsol) + & 
-                            delta(i,k)* (piz_aers_6bands(RH_num(i,k)+1,inu,spsol) - & 
-                            piz_aers_6bands(RH_num(i,k),inu,spsol))
+               piz_ae2b_int = piz_MG_6bands(RH_num(i,k),inu,classbc) +                 &
+                              delta(i,k)* (piz_MG_6bands(RH_num(i,k)+1,inu,classbc) -  &
+                              piz_MG_6bands(RH_num(i,k),inu,classbc))
                    
-             cg_ae2b_int = cg_aers_6bands(RH_num(i,k),inu,spsol) + & 
-                           delta(i,k)* (cg_aers_6bands(RH_num(i,k)+1,inu,spsol) - & 
-                           cg_aers_6bands(RH_num(i,k),inu,spsol))
-
-             tau_ae(i,k,id,inu)    = tmp_var*tau_ae2b_int
-             tau_ae_pi(i,k,id,inu) = tmp_var_pi* tau_ae2b_int
-             piz_ae(i,k,id,inu)    = piz_ae2b_int
-             cg_ae(i,k,id,inu)     = cg_ae2b_int
-                      
+               cg_ae2b_int = cg_MG_6bands(RH_num(i,k),inu,classbc) +                   & 
+                             delta(i,k)* (cg_MG_6bands(RH_num(i,k)+1,inu,classbc) -    &
+                             cg_MG_6bands(RH_num(i,k),inu,classbc))
+
+               tau_ae(i,k,id,inu)    = tmp_var*tau_ae2b_int
+               piz_ae(i,k,id,inu)    = piz_ae2b_int
+               cg_ae(i,k,id,inu)     = cg_ae2b_int
+
+               !--natural aerosols
+               tau_ae2b_int= alpha_MG_6bands(RH_num(i,k),inu,classbc_pi)+                 &
+                             delta(i,k)* (alpha_MG_6bands(RH_num(i,k)+1,inu,classbc_pi) - & 
+                             alpha_MG_6bands(RH_num(i,k),inu,classbc_pi))
+                   
+               piz_ae2b_int = piz_MG_6bands(RH_num(i,k),inu,classbc_pi) +                 &
+                              delta(i,k)* (piz_MG_6bands(RH_num(i,k)+1,inu,classbc_pi) -  &
+                              piz_MG_6bands(RH_num(i,k),inu,classbc_pi))
+                   
+               cg_ae2b_int = cg_MG_6bands(RH_num(i,k),inu,classbc_pi) +                   &
+                             delta(i,k)* (cg_MG_6bands(RH_num(i,k)+1,inu,classbc_pi) -    &
+                             cg_MG_6bands(RH_num(i,k),inu,classbc_pi))
+
+               tau_ae_pi(i,k,id,inu) = tmp_var_pi* tau_ae2b_int
+               piz_ae_pi(i,k,id,inu) = piz_ae2b_int
+               cg_ae_pi(i,k,id,inu)  = cg_ae2b_int
+                        
+             ENDDO
            ENDDO
          ENDDO
-       ENDDO
+
+       !--else treat all other cases of soluble aerosols
+       ELSE
+
+         DO k=1, klev
+           DO i=1, klon
+             tmp_var=m_allaer(i,k,spsol)/1.e6*zdh(i,k)*fac
+             tmp_var_pi=m_allaer_pi(i,k,spsol)/1.e6*zdh(i,k)*fac
+
+             DO inu=1,NSW
+
+               tau_ae2b_int= alpha_aers_6bands(RH_num(i,k),inu,spsol)+ & 
+                             delta(i,k)* (alpha_aers_6bands(RH_num(i,k)+1,inu,spsol) - & 
+                             alpha_aers_6bands(RH_num(i,k),inu,spsol))
+                   
+               piz_ae2b_int = piz_aers_6bands(RH_num(i,k),inu,spsol) + & 
+                            delta(i,k)* (piz_aers_6bands(RH_num(i,k)+1,inu,spsol) - & 
+                              piz_aers_6bands(RH_num(i,k),inu,spsol))
+                   
+               cg_ae2b_int = cg_aers_6bands(RH_num(i,k),inu,spsol) + & 
+                             delta(i,k)* (cg_aers_6bands(RH_num(i,k)+1,inu,spsol) - & 
+                             cg_aers_6bands(RH_num(i,k),inu,spsol))
+
+               tau_ae(i,k,id,inu)    = tmp_var*tau_ae2b_int
+               tau_ae_pi(i,k,id,inu) = tmp_var_pi* tau_ae2b_int
+               piz_ae(i,k,id,inu)    = piz_ae2b_int
+               cg_ae(i,k,id,inu)     = cg_ae2b_int
+                        
+             ENDDO
+           ENDDO
+         ENDDO
+
+         !--external mixture case for soluble BC
+         IF (spsol.EQ.1) THEN
+           piz_ae_pi(:,:,id,:) = piz_ae(:,:,id,:)
+           cg_ae_pi(:,:,id,:)  = cg_ae(:,:,id,:)
+         ENDIF
+
+       ENDIF
         
      ELSE    ! For all aerosol insoluble components
@@ -443,5 +671,5 @@
      DO k=1, klev
        DO i=1, klon
-!--anthropogenic aerosol
+!--all (natural + anthropogenic) aerosol
          tau_allaer(i,k,2,inu)=tau_ae(i,k,id_ASSO4M_phy,inu)+tau_ae(i,k,id_CSSO4M_phy,inu)+ &
                                tau_ae(i,k,id_ASBCM_phy,inu)+tau_ae(i,k,id_AIBCM_phy,inu)+   &
@@ -449,32 +677,35 @@
                                tau_ae(i,k,id_ASSSM_phy,inu)+tau_ae(i,k,id_CSSSM_phy,inu)+   &
                                tau_ae(i,k,id_SSSSM_phy,inu)+ tau_ae(i,k,id_CIDUSTM_phy,inu)
-         tau_allaer(i,k,2,inu)=MAX(tau_allaer(i,k,2,inu),1e-15)
-
-         piz_allaer(i,k,2,inu)=(tau_ae(i,k,id_ASSO4M_phy,inu)*piz_ae(i,k,id_ASSO4M_phy,inu)+ &
-                                tau_ae(i,k,id_CSSO4M_phy,inu)*piz_ae(i,k,id_CSSO4M_phy,inu)+ &
-                                tau_ae(i,k,id_ASBCM_phy,inu)*piz_ae(i,k,id_ASBCM_phy,inu)+ &
-                                tau_ae(i,k,id_AIBCM_phy,inu)*piz_ae(i,k,id_AIBCM_phy,inu)+ &
-                                tau_ae(i,k,id_ASPOMM_phy,inu)*piz_ae(i,k,id_ASPOMM_phy,inu)+ &
-                                tau_ae(i,k,id_AIPOMM_phy,inu)*piz_ae(i,k,id_AIPOMM_phy,inu)+ &
-                                tau_ae(i,k,id_ASSSM_phy,inu)*piz_ae(i,k,id_ASSSM_phy,inu)+ &
-                                tau_ae(i,k,id_CSSSM_phy,inu)*piz_ae(i,k,id_CSSSM_phy,inu)+ &
-                                tau_ae(i,k,id_SSSSM_phy,inu)*piz_ae(i,k,id_SSSSM_phy,inu)+ &
+         tau_allaer(i,k,2,inu)=MAX(tau_allaer(i,k,2,inu),tau_min)
+
+         piz_allaer(i,k,2,inu)=(tau_ae(i,k,id_ASSO4M_phy,inu)*piz_ae(i,k,id_ASSO4M_phy,inu)+   &
+                                tau_ae(i,k,id_CSSO4M_phy,inu)*piz_ae(i,k,id_CSSO4M_phy,inu)+   &
+                                tau_ae(i,k,id_ASBCM_phy,inu)*piz_ae(i,k,id_ASBCM_phy,inu)+     &
+                                tau_ae(i,k,id_AIBCM_phy,inu)*piz_ae(i,k,id_AIBCM_phy,inu)+     &
+                                tau_ae(i,k,id_ASPOMM_phy,inu)*piz_ae(i,k,id_ASPOMM_phy,inu)+   &
+                                tau_ae(i,k,id_AIPOMM_phy,inu)*piz_ae(i,k,id_AIPOMM_phy,inu)+   &
+                                tau_ae(i,k,id_ASSSM_phy,inu)*piz_ae(i,k,id_ASSSM_phy,inu)+     &
+                                tau_ae(i,k,id_CSSSM_phy,inu)*piz_ae(i,k,id_CSSSM_phy,inu)+     &
+                                tau_ae(i,k,id_SSSSM_phy,inu)*piz_ae(i,k,id_SSSSM_phy,inu)+     &
                                 tau_ae(i,k,id_CIDUSTM_phy,inu)*piz_ae(i,k,id_CIDUSTM_phy,inu)) &
                                 /tau_allaer(i,k,2,inu)
-         piz_allaer(i,k,2,inu)=MAX(piz_allaer(i,k,2,inu),0.01)
+         piz_allaer(i,k,2,inu)=MIN(MAX(piz_allaer(i,k,2,inu),0.01),1.0)
+         IF (tau_allaer(i,k,2,inu).LE.tau_min) piz_allaer(i,k,2,inu)=1.0
 
          cg_allaer(i,k,2,inu)=(tau_ae(i,k,id_ASSO4M_phy,inu)*piz_ae(i,k,id_ASSO4M_phy,inu)*cg_ae(i,k,id_ASSO4M_phy,inu)+ &
                                tau_ae(i,k,id_CSSO4M_phy,inu)*piz_ae(i,k,id_CSSO4M_phy,inu)*cg_ae(i,k,id_CSSO4M_phy,inu)+ &
-                               tau_ae(i,k,id_ASBCM_phy,inu)*piz_ae(i,k,id_ASBCM_phy,inu)*cg_ae(i,k,id_ASBCM_phy,inu)+ &
-                               tau_ae(i,k,id_AIBCM_phy,inu)*piz_ae(i,k,id_AIBCM_phy,inu)*cg_ae(i,k,id_AIBCM_phy,inu)+ &
+                               tau_ae(i,k,id_ASBCM_phy,inu)*piz_ae(i,k,id_ASBCM_phy,inu)*cg_ae(i,k,id_ASBCM_phy,inu)+    &
+                               tau_ae(i,k,id_AIBCM_phy,inu)*piz_ae(i,k,id_AIBCM_phy,inu)*cg_ae(i,k,id_AIBCM_phy,inu)+    &
                                tau_ae(i,k,id_ASPOMM_phy,inu)*piz_ae(i,k,id_ASPOMM_phy,inu)*cg_ae(i,k,id_ASPOMM_phy,inu)+ &
                                tau_ae(i,k,id_AIPOMM_phy,inu)*piz_ae(i,k,id_AIPOMM_phy,inu)*cg_ae(i,k,id_AIPOMM_phy,inu)+ &
-                               tau_ae(i,k,id_ASSSM_phy,inu)*piz_ae(i,k,id_ASSSM_phy,inu)*cg_ae(i,k,id_ASSSM_phy,inu)+ &
-                               tau_ae(i,k,id_CSSSM_phy,inu)*piz_ae(i,k,id_CSSSM_phy,inu)*cg_ae(i,k,id_CSSSM_phy,inu)+ &
-                               tau_ae(i,k,id_SSSSM_phy,inu)*piz_ae(i,k,id_SSSSM_phy,inu)*cg_ae(i,k,id_SSSSM_phy,inu)+ &
+                               tau_ae(i,k,id_ASSSM_phy,inu)*piz_ae(i,k,id_ASSSM_phy,inu)*cg_ae(i,k,id_ASSSM_phy,inu)+    &
+                               tau_ae(i,k,id_CSSSM_phy,inu)*piz_ae(i,k,id_CSSSM_phy,inu)*cg_ae(i,k,id_CSSSM_phy,inu)+    &
+                               tau_ae(i,k,id_SSSSM_phy,inu)*piz_ae(i,k,id_SSSSM_phy,inu)*cg_ae(i,k,id_SSSSM_phy,inu)+    &
                                tau_ae(i,k,id_CIDUSTM_phy,inu)*piz_ae(i,k,id_CIDUSTM_phy,inu)*cg_ae(i,k,id_CIDUSTM_phy,inu))/ &
                                (tau_allaer(i,k,2,inu)*piz_allaer(i,k,2,inu))
+         cg_allaer(i,k,2,inu)=MIN(MAX(cg_allaer(i,k,2,inu),0.0),1.0)
 
 !--natural aerosol
+!--ASBCM aerosols take _pi value because of internal mixture option
          tau_allaer(i,k,1,inu)=tau_ae_pi(i,k,id_ASSO4M_phy,inu)+tau_ae_pi(i,k,id_CSSO4M_phy,inu)+ &
                                tau_ae_pi(i,k,id_ASBCM_phy,inu)+tau_ae_pi(i,k,id_AIBCM_phy,inu)+   &
@@ -482,30 +713,32 @@
                                tau_ae_pi(i,k,id_ASSSM_phy,inu)+tau_ae_pi(i,k,id_CSSSM_phy,inu)+   &
                                tau_ae_pi(i,k,id_SSSSM_phy,inu)+ tau_ae_pi(i,k,id_CIDUSTM_phy,inu)
-         tau_allaer(i,k,1,inu)=MAX(tau_allaer(i,k,1,inu),1e-15)
-
-         piz_allaer(i,k,1,inu)=(tau_ae_pi(i,k,id_ASSO4M_phy,inu)*piz_ae(i,k,id_ASSO4M_phy,inu)+ &
-                                tau_ae_pi(i,k,id_CSSO4M_phy,inu)*piz_ae(i,k,id_CSSO4M_phy,inu)+ &
-                                tau_ae_pi(i,k,id_ASBCM_phy,inu)*piz_ae(i,k,id_ASBCM_phy,inu)+ &
-                                tau_ae_pi(i,k,id_AIBCM_phy,inu)*piz_ae(i,k,id_AIBCM_phy,inu)+ &
-                                tau_ae_pi(i,k,id_ASPOMM_phy,inu)*piz_ae(i,k,id_ASPOMM_phy,inu)+ &
-                                tau_ae_pi(i,k,id_AIPOMM_phy,inu)*piz_ae(i,k,id_AIPOMM_phy,inu)+ &
-                                tau_ae_pi(i,k,id_ASSSM_phy,inu)*piz_ae(i,k,id_ASSSM_phy,inu)+ &
-                                tau_ae_pi(i,k,id_CSSSM_phy,inu)*piz_ae(i,k,id_CSSSM_phy,inu)+ &
-                                tau_ae_pi(i,k,id_SSSSM_phy,inu)*piz_ae(i,k,id_SSSSM_phy,inu)+ &
+         tau_allaer(i,k,1,inu)=MAX(tau_allaer(i,k,1,inu),tau_min)
+
+         piz_allaer(i,k,1,inu)=(tau_ae_pi(i,k,id_ASSO4M_phy,inu)*piz_ae(i,k,id_ASSO4M_phy,inu)+   &
+                                tau_ae_pi(i,k,id_CSSO4M_phy,inu)*piz_ae(i,k,id_CSSO4M_phy,inu)+   &
+                                tau_ae_pi(i,k,id_ASBCM_phy,inu)*piz_ae_pi(i,k,id_ASBCM_phy,inu)+  &
+                                tau_ae_pi(i,k,id_AIBCM_phy,inu)*piz_ae(i,k,id_AIBCM_phy,inu)+     &
+                                tau_ae_pi(i,k,id_ASPOMM_phy,inu)*piz_ae(i,k,id_ASPOMM_phy,inu)+   &
+                                tau_ae_pi(i,k,id_AIPOMM_phy,inu)*piz_ae(i,k,id_AIPOMM_phy,inu)+   &
+                                tau_ae_pi(i,k,id_ASSSM_phy,inu)*piz_ae(i,k,id_ASSSM_phy,inu)+     &
+                                tau_ae_pi(i,k,id_CSSSM_phy,inu)*piz_ae(i,k,id_CSSSM_phy,inu)+     &
+                                tau_ae_pi(i,k,id_SSSSM_phy,inu)*piz_ae(i,k,id_SSSSM_phy,inu)+     &
                                 tau_ae_pi(i,k,id_CIDUSTM_phy,inu)*piz_ae(i,k,id_CIDUSTM_phy,inu)) &
                                 /tau_allaer(i,k,1,inu)
-         piz_allaer(i,k,1,inu)=MAX(piz_allaer(i,k,1,inu),0.01)
-
-         cg_allaer(i,k,1,inu)=(tau_ae_pi(i,k,id_ASSO4M_phy,inu)*piz_ae(i,k,id_ASSO4M_phy,inu)*cg_ae(i,k,id_ASSO4M_phy,inu)+ &
-                               tau_ae_pi(i,k,id_CSSO4M_phy,inu)*piz_ae(i,k,id_CSSO4M_phy,inu)*cg_ae(i,k,id_CSSO4M_phy,inu)+ &
-                               tau_ae_pi(i,k,id_ASBCM_phy,inu)*piz_ae(i,k,id_ASBCM_phy,inu)*cg_ae(i,k,id_ASBCM_phy,inu)+ &
-                               tau_ae_pi(i,k,id_AIBCM_phy,inu)*piz_ae(i,k,id_AIBCM_phy,inu)*cg_ae(i,k,id_AIBCM_phy,inu)+ &
-                               tau_ae_pi(i,k,id_ASPOMM_phy,inu)*piz_ae(i,k,id_ASPOMM_phy,inu)*cg_ae(i,k,id_ASPOMM_phy,inu)+ &
-                               tau_ae_pi(i,k,id_AIPOMM_phy,inu)*piz_ae(i,k,id_AIPOMM_phy,inu)*cg_ae(i,k,id_AIPOMM_phy,inu)+ &
-                               tau_ae_pi(i,k,id_ASSSM_phy,inu)*piz_ae(i,k,id_ASSSM_phy,inu)*cg_ae(i,k,id_ASSSM_phy,inu)+ &
-                               tau_ae_pi(i,k,id_CSSSM_phy,inu)*piz_ae(i,k,id_CSSSM_phy,inu)*cg_ae(i,k,id_CSSSM_phy,inu)+ &
-                               tau_ae_pi(i,k,id_SSSSM_phy,inu)*piz_ae(i,k,id_SSSSM_phy,inu)*cg_ae(i,k,id_SSSSM_phy,inu)+ &
+         piz_allaer(i,k,1,inu)=MIN(MAX(piz_allaer(i,k,1,inu),0.01),1.0)
+         IF (tau_allaer(i,k,1,inu).LE.tau_min) piz_allaer(i,k,1,inu)=1.0
+
+         cg_allaer(i,k,1,inu)=(tau_ae_pi(i,k,id_ASSO4M_phy,inu)*piz_ae(i,k,id_ASSO4M_phy,inu)*cg_ae(i,k,id_ASSO4M_phy,inu)+    &
+                               tau_ae_pi(i,k,id_CSSO4M_phy,inu)*piz_ae(i,k,id_CSSO4M_phy,inu)*cg_ae(i,k,id_CSSO4M_phy,inu)+    &
+                               tau_ae_pi(i,k,id_ASBCM_phy,inu)*piz_ae_pi(i,k,id_ASBCM_phy,inu)*cg_ae_pi(i,k,id_ASBCM_phy,inu)+ &
+                               tau_ae_pi(i,k,id_AIBCM_phy,inu)*piz_ae(i,k,id_AIBCM_phy,inu)*cg_ae(i,k,id_AIBCM_phy,inu)+       &
+                               tau_ae_pi(i,k,id_ASPOMM_phy,inu)*piz_ae(i,k,id_ASPOMM_phy,inu)*cg_ae(i,k,id_ASPOMM_phy,inu)+    &
+                               tau_ae_pi(i,k,id_AIPOMM_phy,inu)*piz_ae(i,k,id_AIPOMM_phy,inu)*cg_ae(i,k,id_AIPOMM_phy,inu)+    &
+                               tau_ae_pi(i,k,id_ASSSM_phy,inu)*piz_ae(i,k,id_ASSSM_phy,inu)*cg_ae(i,k,id_ASSSM_phy,inu)+       &
+                               tau_ae_pi(i,k,id_CSSSM_phy,inu)*piz_ae(i,k,id_CSSSM_phy,inu)*cg_ae(i,k,id_CSSSM_phy,inu)+       &
+                               tau_ae_pi(i,k,id_SSSSM_phy,inu)*piz_ae(i,k,id_SSSSM_phy,inu)*cg_ae(i,k,id_SSSSM_phy,inu)+       &
                                tau_ae_pi(i,k,id_CIDUSTM_phy,inu)*piz_ae(i,k,id_CIDUSTM_phy,inu)*cg_ae(i,k,id_CIDUSTM_phy,inu))/ &
                                (tau_allaer(i,k,1,inu)*piz_allaer(i,k,1,inu))
+         cg_allaer(i,k,1,inu)=MIN(MAX(cg_allaer(i,k,1,inu),0.0),1.0)
 
         ENDDO
@@ -513,5 +746,5 @@
     ENDDO
    
-!--waveband 2 and all aerosol
+!--waveband 2 and all aerosol (third index = 2)
   inu=2
   DO i=1, klon
Index: /LMDZ5/trunk/libf/phylmd/rrtm/readaerosol_optic_rrtm.F90
===================================================================
--- /LMDZ5/trunk/libf/phylmd/rrtm/readaerosol_optic_rrtm.F90	(revision 2643)
+++ /LMDZ5/trunk/libf/phylmd/rrtm/readaerosol_optic_rrtm.F90	(revision 2644)
@@ -2,5 +2,5 @@
 !
 SUBROUTINE readaerosol_optic_rrtm(debut, aerosol_couple,  &
-     new_aod, flag_aerosol, itap, rjourvrai, &
+     new_aod, flag_aerosol, flag_bc_internal_mixture, itap, rjourvrai, &
      pdtphys, pplay, paprs, t_seri, rhcl, presnivs, &
      tr_seri, mass_solu_aero, mass_solu_aero_pi, &
@@ -32,4 +32,5 @@
   LOGICAL, INTENT(IN)                      :: new_aod
   INTEGER, INTENT(IN)                      :: flag_aerosol
+  LOGICAL, INTENT(IN)                      :: flag_bc_internal_mixture
   INTEGER, INTENT(IN)                      :: itap
   REAL, INTENT(IN)                         :: rjourvrai
@@ -303,15 +304,17 @@
 !--new aerosol properties
   ! aeropt_6bands for rrtm
-  CALL aeropt_6bands_rrtm( &
-       pdel, m_allaer, pdtphys, rhcl, & 
-       tau_aero, piz_aero, cg_aero,   &
-       m_allaer_pi, flag_aerosol, &
-       zrho ) 
+  CALL aeropt_6bands_rrtm(          &
+       pdel, m_allaer, rhcl,        & 
+       tau_aero, piz_aero, cg_aero, &
+       m_allaer_pi, flag_aerosol,   &
+       flag_bc_internal_mixture, zrho ) 
 
   ! aeropt_5wv only for validation and diagnostics
-  CALL aeropt_5wv_rrtm(                    &
-       pdel, m_allaer,                &
-       pdtphys, rhcl, aerindex,       & 
-       flag_aerosol, pplay, t_seri,   &
+  CALL aeropt_5wv_rrtm(              &
+       pdel, m_allaer,               &
+       rhcl, aerindex,               & 
+       flag_aerosol,                 & 
+       flag_bc_internal_mixture,     & 
+       pplay, t_seri,                &
        tausum_aero, tau3d_aero )
 
