In [1]:
import numpy as np
import escape as esc

esc.require("0.9.8")
Loading material database from C:\dev\escape-core\python\src\escape\scattering\..\data\mdb\materials.db

Functors¶

The intensity equation in scattering depends on sample type (thin film, bulk, isotropic or anisotropic), setup (detector type, slits, beam profile), and often diffuse scattering or background terms that are included as polynomials or other expressions. In real experiments, these parts may represent absorption, resolution broadening, baseline signal, or the response of the instrument.

To support custom models and laboratory-specific setups in one framework, ESCAPE uses a flexible functor object: an object that behaves like a function of one or more variables (up to five). Functors support arithmetic and standard mathematical functions so that model formulas stay close to algebraic notation. They can be evaluated at single points or over coordinate arrays, which makes them useful both for plotting and for fitting experimental data.

Key building blocks are variables (the arguments of the function), parameters (optimizable or fixed constants), and constant functors (functions that always return a given value). Combined with operators and math functions, they allow arbitrary expressions. Below we introduce variables, constant functors, and expressions.

Variables¶

Variables are created with esc.var; they represent the arguments of a functor. For example:

In [2]:
X = esc.var("X")
Y = esc.var("Y")
Z = esc.var("Z")

Variables support algebraic expressions and operators. If one writes $x+y$, the result of this expression is a two-dimensional functor, i.e. function with two variables $f(x, y)$

In [3]:
f = X + Y
print("Type: ", type(f))
print(f.num_variables)
print(f.variables)
Type:  <class 'escape.core.entities.functor_obj'>
2
[variable(name='X'), variable(name='Y')]

Const functor¶

A constant functor always returns a given value (a number or a parameter), regardless of the variable values at which it is evaluated. You construct it with esc.func: the first argument is a variable or list of variables, the second is the constant (float or parameter). For one variable use a single variable; for several variables pass a list, e.g. esc.func([X, Y], c).

In [4]:
p1 = esc.par("par1", 23.5)

# one-dimensional const functors
f1 = esc.func(X, p1)
f2 = esc.func(Y, 10)

# two-dimensional const functors
f3 = esc.func([X, Y], p1)

print(f1(10))
print(f2(10))
print(f3(1, 2))
23.5
10.0
23.5

Expressions¶

Functors support arithmetic operators and mathematical functions: sin(), cos(), tan(), sinh(), cosh(), tanh(), exp(), sqrt(), log(), log10(), pow(), min(), max(), erf(). You can combine variables and parameters into expressions; the result is a functor whose number of variables is determined by the variables that appear in it. For use in optimizers, esc.reduce can bind some variables to parameters or constants to reduce the number of variables. Below we show a simple oscillating functor with damping amplitude.

In [5]:
q = esc.par("Damping", 0.05, userlim=[0, 2])
expr = esc.pow(esc.sin(X), 2.0) * esc.exp(-q * esc.abs(X))

expr.show().config(title="Damping sine")
Out[5]:
In [6]:
# two-dimensional example
Qx = esc.par("Damping X", 0.05, userlim=[0, 2])
Qy = esc.par("Damping Y", 0.05, userlim=[0, 2])

expr2d = (
    esc.pow(esc.cos(X), 2.0)
    * esc.pow(esc.cos(Y), 2.0)
    * esc.exp(-Qx * esc.abs(X))
    * esc.exp(-Qy * esc.abs(Y))
)

x = np.linspace(-5, 5, 100)
y = np.linspace(-5, 5, 100)

# create meshgrid of both coordinates
xv, yv = np.meshgrid(x, y)
# create coordinates array of type [x0, y0, x1, y1,...]
coords = np.column_stack([xv.flatten(), yv.flatten()])

# for the 2d case one has to provide a plot_type, coordinates index and number of rows
# cbmin, cbmax and cblof are the colorbox properties

expr2d.show(coordinates=coords).config(title="Map Plot")
Out[6]:

Evaluation: scalar and array¶

A functor can be called in two ways. Scalar call: pass one value per variable (e.g. f(0.5) for one variable, f(0.5, 0.3) for two); the result is a single float. Array call: pass a 1D NumPy array (or double_array_obj) of coordinates; the functor is evaluated at each point in parallel and returns a new float64 array, or you can pass a second array to fill in-place: f(x_array, y_array). For multi-variable functors, the input array must contain stacked coordinates (e.g. for 2D, length is n×2: x0,y0, x1,y1, ...). Array evaluation is used internally when fitting (model cost over data points) and can be controlled with num_threads and device (see below).

In [7]:
# Scalar vs array evaluation
F = esc.sin(X) * esc.exp(-0.1 * X)
print("Scalar F(1.0) =", F(1.0))

x_arr = np.linspace(0, 5, 100)
y_arr = F(x_arr)  # returns new array
print("Array length:", len(y_arr))

y_out = np.empty(100)
F(x_arr, y_out)  # in-place; y_out is filled
print("In-place matches:", np.allclose(y_arr, y_out))

# Thread count (0 = default/auto); device is "cpu" or a ZeroMQ endpoint for remote evaluation
F.num_threads = 4
print("num_threads:", F.num_threads, "device:", F.device)
Scalar F(1.0) = 0.7613944332457532
Array length: 100
In-place matches: True
num_threads: 4 device: cpu
C:\Users\User\AppData\Local\Temp\ipykernel_34128\1578388116.py:6: UserWarning:

Input array does not own its data (e.g. it is a view or slice); data will be copied

C:\Users\User\AppData\Local\Temp\ipykernel_34128\1578388116.py:10: UserWarning:

Input array does not own its data (e.g. it is a view or slice); data will be copied

Scale functor¶

esc.scale(func, x=(x1, x2), y=(y1, y2)) returns a new functor that is an affine transformation of the one-variable functor func so that it passes through the points (x1, y1) and (x2, y2). Useful for normalizing a curve to a given range or for fitting with fixed end points.

In [8]:
# Scale a functor to pass through (0, 0) and (1, 100)
base = esc.sin(X)
scaled = esc.scale(base, x=(0, np.pi / 2), y=(0, 100))
print("base(0) =", base(0), ", scaled(0) =", scaled(0))
print("base(pi/2) =", base(np.pi / 2), ", scaled(pi/2) =", scaled(np.pi / 2))
base(0) = 0.0 , scaled(0) = 0.0
base(pi/2) = 1.0 , scaled(pi/2) = 100.0

Variables order¶

ESCAPE supports expressions with functors of different domain size and preserves the order of variables in the resulted functor. For example, $f_1(z, x)+f_2(x)+f_3(y)$ will result in the following functor $f_r(z, x, y)$, not $f_r(x, y, z)$ as maybe desired. One has to keep this in mind when creating large complex expressions. One can always check the order of variables using variables property of functors.

In [9]:
F1 = esc.func([Z, X], 10, name="F1(z,x)")
F2 = esc.func(X, 0.1, name="F2(x)")
F3 = esc.func(Y, 0.05, name="F3(y)")

Fr = F1 + F2 + F3

print(Fr.variables)
[variable(name='X'), variable(name='Y'), variable(name='Z')]

Boolean functors and conditional expressions¶

If two functors are compared using one of the standard comparison operators >, <, >=, <=, ==, !=, a product of this comparison is a boolean functor. For example, we have two functors f1(x)=x^2 and f2(x,y)=x+y

In [10]:
F1 = X * X
F2 = X + Y
In [11]:
Fb = F1 >= F2
In [12]:
print(type(Fb))
Fb.variables
<class 'escape.core.entities.bool_functor_obj'>
Out[12]:
[variable(name='X'), variable(name='Y')]

Fb is a boolean functor of two variables $f(x, y)$. When beeing called with with values of variables $x$ and $y$, it calls f1(x) and f2(x, y) and returns a result of f1(x) >= f2(x). Boolean functors support logical operators: or - |, and - & and not - ~.

In [13]:
Fb(1, 10)
Out[13]:
False
In [14]:
Fb(4, 5)
Out[14]:
True

Boolean functors are typically used with conditional expressions. The first argument of the conditional method is a condition itself, i.e. a boolean functor, the second one is a functor which has to be executed if the condition returns true, and the last one is executed if the condition returns false.

In [15]:
P = esc.par("P", 10, userlim=[-1e-3, 1e-3])
B = esc.par("B", 0, userlim=[-10, 10])
F1 = B * X + P
F2 = P * X + B
Fc = esc.conditional((X >= -5) & (X <= 5), F1, F2)
In [16]:
a = np.arange(-10, 10, 0.01)
p1 = F1.show().config(title="F1(x)")
p2 = F2.show().config(title="F2(x)")
p3 = Fc.show().config(title="Fc(x)")
esc.show(p1, p2, p3).config(title="Conditional expression")
Out[16]:

Functors with complex return type¶

In scattering the connection between experiment geometry and equations describing measured intensity is done in terms of reciprocal space. In the first approximation, called Born approximation, the conversion between real- and reciprocal- space is performed using Fourier Transform. The result of Fourier Transform can be a complex function. Thus, it make sense to introduce a functor with complex return type.

In [17]:
# Constant complex functor
Cf = esc.cfunc([X], -1j)
In [18]:
type(Cf)
Out[18]:
escape.core.entities.cplx_functor_obj

Complex functors support the same arithmetic operators and mathematical functors as functors with double return type. Additionally there are several other methods like, norm, real, imag, conjugate which require cplx_functor_obj input parameter. Few examples are given below.

In [19]:
F = esc.real((esc.exp(1j * X) + esc.exp(-1j * X)) / 2.0)
In [20]:
F.show()
Out[20]:
In [21]:
F = esc.real((esc.exp(1j * X) - esc.exp(-1j * X)) / (2 * 1j))
In [22]:
F.show()
Out[22]:
In [23]:
F = esc.imag(esc.sin(-1j * X)) + esc.real(esc.cos(-1j * X))
In [24]:
F.show()
Out[24]:

Functor minimizer¶

There is a special type of functor used to find the minimum of another functor. When called, it minimizes the functor parameters for the given variable values.

For example:

your task is to find a minimum of the following functor:

$f= sin(3x)+sin(1/3x)$

In [25]:
# Let's create a functor

a = esc.par("a", 1, userlim=[0, 2 * np.pi])
f0 = esc.sin(3 * X) + esc.sin(X * 1 / 3)
display(f0.show(coordinates=np.linspace(0, 2 * np.pi, 100)))
# since the function minimizer finds a minimum of the function
# by optimizing its parameters, we need to substitute the variable
# with the parameter a. This is done using the reduce method.
f = esc.reduce(f0, X, a)

# Now we create a minimizer using differential evolution algorithm
fmin = esc.diffevol(f, ftol=1e-6, maxiter=1000, popsize=15, status_exc=True)

print("f=", f(), "a=", a.value)
print("fmin=", fmin(), "amin=", a.value)
f= 0.46831470485601945 a= 1.0
fmin= -0.504661876460889 amin= 1.5384760658368035
In [26]:
# Let's add a constraint to the parameter a
f.constrain((a > 3) & (a <= 5.5))

# Now we can find a minimum of the function with constraint
print("fmin=", fmin(), "amin=", a.value)
fmin= -0.06103802623550969 amin= 3.652371976970784
In [27]:
# Now for 2d functor
a = esc.par("a", 1, userlim=[0, 2 * np.pi])
b = esc.par("b", 1, userlim=[0, 2 * np.pi])
f0 = esc.sin(3 * X) + esc.sin(X * 1 / 3) + esc.sin(3 * Y) + esc.sin(Y * 1 / 3)
f = esc.reduce(f0, X, a, Y, b)

xc = np.linspace(0, 2 * np.pi, 100)
yc = np.linspace(0, 2 * np.pi, 100)

xv, yv = np.meshgrid(xc, yc)
coords = np.column_stack([xv.flatten(), yv.flatten()])

f0.show(coordinates=coords)
Out[27]:
In [28]:
# Create a minimizer and call it
fmin = esc.diffevol(f, ftol=1e-6, maxiter=1000, popsize=15, status_exc=True)

print("f=", f(), "a=", a.value, "b=", b.value)
print("fmin=", fmin(), "amin=", a.value, "bmin=", b.value)
f= 0.9366294097120389 a= 1.0 b= 1.0
fmin= -1.0093237521435645 amin= 1.5384698170193163 bmin= 1.5384594843272343
In [29]:
# Add a constraint to the parameters
f.constrain(a**2 + b**2 >= 5)

print("fmin=", fmin(), "amin=", a.value, "bmin=", b.value)
fmin= -0.9930720633587865 amin= 1.5811793289463358 bmin= 1.581098331845101
In [ ]:
 
In [ ]: