Implement deterministic nontrivial function of random variables

Questions
5

Solved! When I was preparing the submission I found the error in the matrix dimensions.

The complete code is here,

import numpy as np
import pymc3 as pm
import theano.tensor as tt

with pm.Model() as model:
        a = pm.Uniform('a', 0, 60)
        R = pm.Uniform('R', 0, 60)
        N0 = pm.Uniform('N0', -2, 2)
        f = pm.Uniform('f', -1, 1)
        d = pm.Uniform('d', 0, 50)
        B0 = pm.Uniform('B0', 0, 1000)
        alpha = pm.Uniform('alpha', -np.pi/2, np.pi/2)
        xc = pm.Uniform('xc', 0, 113)
        yc = pm.Uniform('yc', 0, 53)

sz1 = 112 #113
sz2 = 52 #53
nx, ny = (sz1, sz2)
x1 = np.linspace(0, sz1, sz1+1)
y1 = np.linspace(0, sz2, sz2+1)
xv, yv = np.meshgrid(x1, y1)

def torusmodel(inputParam):
    a, R, N0, f, d, B0, alpha, xc, yc = inputParam
    
    x=(xv-xc+0.5)*tt.cos(alpha) + (yv-yc+0.5)*tt.sin(alpha)
    y=-(xv-xc+0.5)*tt.sin(alpha) + (yv-yc+0.5)*tt.cos(alpha)
    xr=tt.sqrt(x**2+d**2) - R
    rho=tt.sqrt(xr**2 + y**2)
    u=tt.sqrt(x**2 + d**2)
    costh=xr/rho
    Nt=N0*(1+f*(rho/a)**2)
    mag = (x*(1+g*(costh)) - 2*Nt*d*y/u)*(-1)*B0*tt.exp((-1)*(rho/a)**2)/u

    return mag

inputParam = (10,30,0.5,1,0,15,1000,10,50,30)

with model:
    torusmodel_pm = pm.Deterministic('torusmodel_pm', torusmodel(inputParam))

with model:
    observations = pm.Normal( "obs",  mu=torusmodel_pm, tau=1/0.1, observed=data.values)

with model:    
    trace = pm.sample(10000)

Auto-assigning NUTS sampler...
Initializing NUTS using jitter+adapt_diag...
Multiprocess sampling (4 chains in 4 jobs)
NUTS: [yc, xc, alpha, B0, d, f, N0, R, a]
Sampling 4 chains: 100%|██████████| 42000/42000 [05:55<00:00, 118.21draws/s]

in which data.values is a (53, 113) real matrix of solar magnetic field observations. It compiles now, although the traces for all the parameters doesn’t seem to converge (see example for parameter ‘R’ below).
R_trace

But this is of course another problem.

Thank you for the feedback!
Francisco