Skip to content

Commit

Permalink
avoid "is 0" and "is not 0"
Browse files Browse the repository at this point in the history
  • Loading branch information
ghislainp committed Jun 27, 2023
1 parent 314322b commit 4796e6f
Show file tree
Hide file tree
Showing 2 changed files with 22 additions and 18 deletions.
14 changes: 9 additions & 5 deletions smrt/core/lib.py
Original file line number Diff line number Diff line change
Expand Up @@ -168,7 +168,7 @@ class smrt_matrix(object):

def __init__(self, mat, mtype=None):

if mat is 0:
if is_zero_scalar(mat):
self.values = np.float64(0.) # 0, but can be used as a numpy thing
self.mtype = "0"
else:
Expand Down Expand Up @@ -217,8 +217,8 @@ def full(dims, value, mtype=None):
def npol(self):
return self.values.shape[0]

def isnull(self):
return isnull(self)
def is_equal_zero(self):
return is_equal_zero(self)

def compress(self, mode=None, auto_reduce_npol=False):
"""compress a matrix. This comprises several actions:
Expand Down Expand Up @@ -348,13 +348,17 @@ def __repr__(self):
return str("smrt_matrix %s %s" % (self.mtype, shape)) + "\n" + str(self.values)


def isnull(m):
def is_zero_scalar(m):
return np.isscalar(m) and (m == 0)


def is_equal_zero(m):
"""return true if the smrt matrix is null"""

if isinstance(m, smrt_diag):
m = m.diagonal()

return (m is 0) or \
return is_zero_scalar(m) or \
(getattr(m, "mtype", None) == "0") or \
(not np.any(m))

Expand Down
26 changes: 13 additions & 13 deletions smrt/rtsolver/dort.py
Original file line number Diff line number Diff line change
Expand Up @@ -27,7 +27,7 @@
# local import
from ..core.error import SMRTError, smrt_warn
from ..core.result import make_result
from smrt.core.lib import smrt_matrix, smrt_diag, isnull
from smrt.core.lib import smrt_matrix, smrt_diag, is_equal_zero
from smrt.core.optional_numba import numba
# Lazy import: from smrt.interface.coherent_flat import process_coherent_layers

Expand Down Expand Up @@ -429,25 +429,25 @@ def dort_modem_banded(self, m, streams, eigenvalue_solver,
Tbottom_lp1 = interfaces.transmission_bottom(l, m, compute_coherent_only)
# the size of Tbottom_lp1 can be the nsl_npol in general or nslp1_npol if only the specular is present
# and some streams are subject to total reflection.
if not isnull(Tbottom_lp1):
if not is_equal_zero(Tbottom_lp1):
ns_npol_common_bottom = min(Tbottom_lp1.shape[0], nslp1_npol)
todiag(bBC, il_top[l + 1], j, -matmul(Tbottom_lp1, Ed, transb)[:ns_npol_common_bottom, :])

# fill the vector
if m == 0 and self.temperature is not None and self.temperature[l] > 0:
if isnull(Rtop_l):
if is_equal_zero(Rtop_l):
b[il_topl:il_topl + nsl_npol, :] -= self.temperature[l] # to be put at layer (l)
else:
b[il_topl:il_topl + nsl_npol, :] -= ((1.0 - muleye(Rtop_l)) * self.temperature[l])[:, np.newaxis] # a mettre en (l)
# the muleye comes from the isotropic emission of the black body

if l < nlayer - 1 and self.temperature[l] > 0 and not isnull(Tbottom_lp1):
if l < nlayer - 1 and self.temperature[l] > 0 and not is_equal_zero(Tbottom_lp1):
b[il_top[l + 1]:il_top[l + 1] + ns_npol_common_bottom, :] += \
(muleye(Tbottom_lp1) * self.temperature[l])[:ns_npol_common_bottom, np.newaxis] # to be put at layer (l + 1)

if l == 0: # Air-snow interface
Tbottom_air_down = interfaces.transmission_bottom(-1, m, compute_coherent_only)
if not isnull(Tbottom_air_down):
if not is_equal_zero(Tbottom_air_down):
ns_npol_common_bottom = min(Tbottom_air_down.shape[0], nsl_npol) # see the comment on Tbottom_lp1
b[il_topl:il_topl + ns_npol_common_bottom, :] += matmul(Tbottom_air_down, intensity_down_m)

Expand All @@ -462,26 +462,26 @@ def dort_modem_banded(self, m, streams, eigenvalue_solver,

if l > 0:
Ttop_lm1 = interfaces.transmission_top(l, m, compute_coherent_only)
if not isnull(Ttop_lm1):
if not is_equal_zero(Ttop_lm1):
ns_npol_common_top = min(Ttop_lm1.shape[0], nslm1_npol) # see the comment on Tbottom_lp1
todiag(bBC, il_bottom[l - 1], j, -matmul(Ttop_lm1, Eu, transt)[:ns_npol_common_top, :]) # to be put at layer (l - 1)

# fill the vector
if m == 0 and self.temperature is not None and self.temperature[l] > 0:
if isnull(Rbottom_l):
if is_equal_zero(Rbottom_l):
b[il_bottoml:il_bottoml + nsl_npol, :] -= self.temperature[l] # to be put at layer (l)
else:
b[il_bottoml:il_bottoml + nsl_npol, :] -= \
((1.0 - muleye(Rbottom_l)) * self.temperature[l])[:, np.newaxis] # to be put at layer (l)
if l > 0 and not isnull(Ttop_lm1):
if l > 0 and not is_equal_zero(Ttop_lm1):
b[il_bottom[l - 1]:il_bottom[l - 1] + ns_npol_common_top, :] += \
(muleye(Ttop_lm1) * self.temperature[l])[:ns_npol_common_top, np.newaxis] # to be put at layer (l - 1)

if m == 0 and l == nlayer - 1 and self.snowpack.substrate is not None and \
self.snowpack.substrate.temperature is not None and self.temperature is not None:
Tbottom_sub = interfaces.transmission_bottom(l, m, compute_coherent_only)
ns_npol_common_bottom = min(Tbottom_sub.shape[0], nsl_npol) # see the comment on Tbottom_lp1
if not isnull(Tbottom_sub):
if not is_equal_zero(Tbottom_sub):
b[il_bottoml:il_bottoml + ns_npol_common_bottom, :] += \
(muleye(Tbottom_sub) * self.snowpack.substrate.temperature)[:ns_npol_common_bottom, np.newaxis] # to be put at layer (l)

Expand Down Expand Up @@ -533,7 +533,7 @@ def muleye(x):

if isinstance(x, smrt_diag):
return x.diagonal()
elif (x is 0) or (len(x.shape) == 0):
elif (x.is_zero_scalar()) or (len(x.shape) == 0):
return np.atleast_1d(x)
else:
assert len(x.shape) == 2
Expand Down Expand Up @@ -645,7 +645,7 @@ def solve(self, m, compute_coherent_only, debug_A=False):
# calculate the A matrix. Eq (12), or 0 if compute_coherent_only
A = self.ft_even_phase.compress(mode=m, auto_reduce_npol=True) if not compute_coherent_only else 0

if isnull(A):
if is_equal_zero(A):
# the solution is trivial
beta = invmu * self.ke(mu, npol=npol).compress().diagonal()
E = np.eye(2 * n, 2 * n)
Expand Down Expand Up @@ -898,7 +898,7 @@ def combine_coherent_diffuse_matrix(coh, diff, m, compute_coherent_only):

mat_coh = coh.compress(mode=m, auto_reduce_npol=True)

if (not compute_coherent_only) and (diff is not 0) and (not diff.isnull()):
if (not compute_coherent_only) and (not is_equal_zero(diff)):
# the coef comes from the integration of \int dphi' cos(m (phi-phi')) cos(n phi')
# m=n=0 --> 2*np.pi
# m=n > 1 --> np.pi
Expand All @@ -916,7 +916,7 @@ def combine_coherent_diffuse_matrix(coh, diff, m, compute_coherent_only):


def normalize_diffuse_matrix(mat, mu_st, mu_i, weights):
if mat.isnull():
if is_equal_zero(mat):
return mat

if mat.mtype == "dense5":
Expand Down

0 comments on commit 4796e6f

Please sign in to comment.