coordinate dependent forcing for 3D Navier Stokes

[Makrand Ajay Khanwale]

Hello developers and users,

I am trying to implement a coordinate dependent linear forcing to 3D Navier Stokes.

It looks like `A_ij(y)u_j`  where A_ij is a matrix u_j is velocity. The A_ij are zero in the middle of the domain and smoothly go to their values near the wall (y-direction).  Therefore, an if else condition on the `y` coordinate and looping over coordinates is required. 

I am trying to do this in the General Function framework but I keep running into errors (I am guessing I am not correctly using the API of dedalus).  Pasting my general function setup just as a reference.  I would really appreciate suggestions on better ways to do this or corrections to the way I am using the API. Really appreciate the help.  The function below has no `y` dependence currently, but I would like to add that.

```# Define forcing
f = dist.Field(bases=(xbasis, ybasis, zbasis))
domain = f.domain
tensorsig_f = f.tensorsig
def forcing_grid_function(*args):
    y=args[0]["g"]
    #u_in["g"][0]=args[1]["g"]
    #u_in["g"][1]=args[2]["g"]
    #u_in["g"][2]=args[3]["g"]
    #u_hp_filtered = u_in
    args[1].high_pass_filter(0.50)
    args[2].high_pass_filter(0.50)
    args[3].high_pass_filter(0.50)
   
    Ymax = np.max(y)
    Ymin = np.min(y)
    nearWallLengthScale = Ly/8
    a_force = np.array([1.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 1.0])
    forcing_u['g'][0] = F * (a_force[0]*args[1] + a_force[1]*args[2] + a_force[2]*args[3])
    forcing_u['g'][1] = F * (a_force[3]*args[1] + a_force[4]*args[2] + a_force[5]*args[3])
    forcing_u['g'][2] = F * (a_force[6]*args[1] + a_force[7]*args[2] + a_force[8]*args[3])

    print(forcing_u)
    return forcing_u


def forcing_u(*args, domain=domain, F=forcing_grid_function):
    return d3.GeneralFunction(
    dist=dist,
    domain=domain,
    tensorsig=(coords,),
    dtype=np.float64,
    layout="g",
    func=F,
    args=args,
)
```

[Miguel Beneitez]

Hi Makrand,

is this being done in serial? I think that if you are trying to do this in parallel probably you want to use the local y in each processor (see manipulating field data here: https://dedalus-project.readthedocs.io/en/latest/notebooks/dedalus_tutorial_2.html#2.1:-Fields). 

Also if your matrix A does not change in time, you can define this before defining your problem and it directly there without needing to use a general function (in my opinion). 

It would also be helpful if you can provide a minimal example of this and what kind of errors you might be getting. 

All the best,

Miguel

--
You received this message because you are subscribed to the Google Groups "Dedalus Users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to dedalus-user...@googlegroups.com.
To view this discussion visit https://groups.google.com/d/msgid/dedalus-users/ddac9f8e-0be6-4ffd-a588-5ec63ab4419en%40googlegroups.com.