simpeg / simpeg

Simulation and Parameter Estimation in Geophysics - A python package for simulation and gradient based parameter estimation in the context of geophysical applications.
http://simpeg.xyz
MIT License
496 stars 262 forks source link

ENH: physical properties as functions? #1345

Open jcapriot opened 7 months ago

jcapriot commented 7 months ago

Proposed new feature or change:

How should we handle physical properties that are themselves functions of parameters that are not the model (aka a physical property that is dependent on how it is observed).

As a concrete question: How should we account for time (or frequency) dependent conductivity?

Clearly a simple solution for the frequency dependent case is to have a single simulation for each frequency we are interested in, (requiring much duplication of internally cached items) but I'm interested in a more general solution.

Say our conductivity is a function of time (as in there might be some induced polarization effects) and we want to simulate time-domain EM data. There is not a good way for us to do that at the moment directly (it would involve creating subclasses of the TDEM simulations, each with their own respective form of the time-dependent conductivity, similar to how the SIP class is specifically using a stretched exponential to handle the time-dependent conductivity).

What if instead of the physical properties being arrays (even in the case of defining a physical property map they evaluate to be arrays), we also allowed them to be functions of some tertiary variable that would then return the correct sigma array.

i.e.

def sigma_func(omega):
    return sigma_0 * some_fancy_cole_cole_model(omega)

sim = SimulationTime(sigma=sigma_func)

Also somehow combine it with a map:

def sigma_map_func(m, time):
    return maps.Exp(m) * some_fancy_cole_cole_model(omega)

Maybe it might be better to elevate this to a dedicated class type for more reliable checking if the physical property is one of these proposed things:

class PropertyFunction:
    def __init__(self, function, ...)  # not quite sure what all would go here
        pass

    # might also need to define some derivatives of the function for certain simulations?
    # e.g. where previous time-independence was assumed to move conductivity outside a time derivative.

Function dependent mapping would likely require some sort of dedicated mapping class:

class FunctionMapping:

    def __init__(self, ...):

    def deriv_m(x, m):
        pass
    ...

# likely used as
sigma = mapping_func(omega, model)
sigma_m_deriv = mapping_func.deriv_m(omega, model)
sigma_t_deriv = mapping_func.deriv_x(omega, model)  # Is this necessary? (Probably)

I realize that this is probably a large project involving a lot of changes and selection cases. It would go into a lot of the fundamental mass matrix code and likely require some reformulations for previous simulations that assumed time-independence of a property, but that shouldn't necessarily be a reason to ignore this forever!

My intention would be to use this page as a springboard to discuss potential ideas details and roadblocks.

ghwilliams commented 7 months ago

It makes sense. It looks very natural to treat the physical properties in the same foot as other quantities like the physical fields. Maybe, to not create another hierarchy of base classes, physical properties could inherit from Fields related classes. Just an idea, if physical properties would be treated as just another type of fields evaluated in the simulation grid.