Diffusion with an "Arrhenius" source term

[Jason Brown]

Hi,

I have a 1D diffusion problem with a source term that is exponential in two of the variables.  It's modeling a chemical reaction with an Arrhenius-type rate relation - see attached image.

When I try to simulate this - see attached working example - I get the following error:

TypeError: Cannot cast array data from dtype('O') to dtype('float64') according to the rule 'same_kind'

Full traceback:

Traceback (most recent call last):

  File "/Users/jasonbrown/Documents/02Work/GDOT Mass Concrete/CodesDedalus/GF_minimum_example.py", line 153, in <module>

    solver.step(dt)

  File "/Users/jasonbrown/opt/anaconda3/envs/dedalus/lib/python3.8/site-packages/dedalus/core/solvers.py", line 507, in step

    self.timestepper.step(self, dt)

  File "/Users/jasonbrown/opt/anaconda3/envs/dedalus/lib/python3.8/site-packages/dedalus/core/timesteppers.py", line 570, in step

    evaluator.evaluate_scheduled(**evaluator_kw)

  File "/Users/jasonbrown/opt/anaconda3/envs/dedalus/lib/python3.8/site-packages/dedalus/core/evaluator.py", line 107, in evaluate_scheduled

    self.evaluate_handlers(scheduled_handlers, wall_time=wall_time, sim_time=sim_time, iteration=iteration, **kw)

  File "/Users/jasonbrown/opt/anaconda3/envs/dedalus/lib/python3.8/site-packages/dedalus/core/evaluator.py", line 142, in evaluate_handlers

    tasks = self.attempt_tasks(tasks, id=id)

  File "/Users/jasonbrown/opt/anaconda3/envs/dedalus/lib/python3.8/site-packages/dedalus/core/evaluator.py", line 187, in attempt_tasks

    output = task['operator'].attempt(**kw)

  File "/Users/jasonbrown/opt/anaconda3/envs/dedalus/lib/python3.8/site-packages/dedalus/core/future.py", line 175, in attempt

    return self.evaluate(id=id, force=False)

  File "/Users/jasonbrown/opt/anaconda3/envs/dedalus/lib/python3.8/site-packages/dedalus/core/future.py", line 132, in evaluate

    a_eval = a.evaluate(id=id, force=force)

  File "/Users/jasonbrown/opt/anaconda3/envs/dedalus/lib/python3.8/site-packages/dedalus/core/future.py", line 161, in evaluate

    self.operate(out)

  File

"/Users/jasonbrown/opt/anaconda3/envs/dedalus/lib/python3.8/site-packages/dedalus/core/operators.py", line 227, in operate

    np.copyto(out.data, self.func(*self.args, **self.kw))

  File "<__array_function__ internals>", line 5, in copyto

TypeError: Cannot cast array data from dtype('O') to dtype('float64') according to the rule 'same_kind'

I've seen a similar error in other messages here but I haven't been able to decipher enough to fix the problem. (Side note: I'm using a GeneralFunction to implement taking the exponential)

Any help is appreciated, and thanks,

Jason

for the model).

[Geoffrey Vasil]

Hi Jason, 

If you PDF equations match you Dedalus code, then I think I see an issue. 

You need to move the nonlinear terms to the right-hand side. Dedalus time steps left-side terms implicitly, and right terms explicitly. 

But the implicit solve only works for linear terms. 

Hope this helps. 

-Geoff

Sent from my iPhone

On 19 May 2021, at 6:09 am, Jason Brown <jbri...@gmail.com> wrote:

Hi,

<Basic model.png>

for the model).

--
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 on the web visit https://groups.google.com/d/msgid/dedalus-users/7f7e5567-1ae0-4e0f-91e8-6b05857fe164n%40googlegroups.com.

<GF_minimum_example.py>

<Basic model.png>

[Keaton Burns]

Hi Jason,

When you pass “e” as a problem parameter, it gets internally casted to a dedalus Scalar object.  It’s actual numerical value should then be accessed inside the general function using the “.value” attribute, i.e. the “f_power” function should be modified to:

def f_power(*args):

    base = args[0].value

    exponent = args[1]['g']

    return np.power(base, exponent)

It looks like this prevents the TypeError, but the script goes unstable after a few iterations.

-Keaton

To view this discussion on the web visit https://groups.google.com/d/msgid/dedalus-users/8710E6FC-89BD-4032-8913-04C22B105150%40gmail.com.

[Jason Brown]

Thanks Geoff, good catch - it was my PDF that was wrong and the code has the nonlinear terms on the RHS.

To view this discussion on the web visit https://groups.google.com/d/msgid/dedalus-users/8710E6FC-89BD-4032-8913-04C22B105150%40gmail.com.

[Jason Brown]

Thanks Keaton, this fixed the problem.

My ‘real’ script (non minimal example) works and is stable - my minimum example may (?) go haywire since I stripped out some material properties to make it simpler to read, but didn’t scale the rest of the minimum example model…(??).

Jason

To view this discussion on the web visit https://groups.google.com/d/msgid/dedalus-users/CADZXxBj_scg1cYE2cSyQsHsW-z7ogLsuXefXH_tXdXDYnB1KqA%40mail.gmail.com.