MATLAB Can COMSOL Integrate MATLAB Functions for Reaction-Diffusion Models?

[TWNien]

Hello,

Now I am working on the convection and diffusion problem with chemical reactions. The question is simple (but the answer is hard to me...) : I want to put a MATLAB function in the generation term, ex: the original one is " -k*c" , and I write a function r_c(x) which returns -k*x. In this case c is the dependent variable, and k is the reaction constant.

It seems to me that COMSOL can't do the Jacobian assembly for m-function with dependent variable (don't know how to compute the derivatives of m-function).

Does anyone have the same experience with the idea? Thanks for sharing.

Daniel

 

[PerennialII]

Yeah, if understood about right that's a problem haven't gotten around to myself (when using m - file - specified functions within the equation systems) directly but rather by specifying the derivatives manually (in the function definition). Often by specifying other m - file functions, where in complex cases compute the derivatives for example using numerical differentiation (in case analytical doesn't exist).

 

[TWNien]

Yeah, if understood about right that's a problem haven't gotten around to myself (when using m - file - specified functions within the equation systems) directly but rather by specifying the derivatives manually (in the function definition). Often by specifying other m - file functions, where in complex cases compute the derivatives for example using numerical differentiation (in case analytical doesn't exist).

Wow, thanks for your QUICK reply. I am so impressed.

Do you have any idea that other software can do this? I've tried specified other m-file for the derivatives, but not succeed yet even with the simplest problem. woe... There is still a long way to go if I want to incorporate more reaction rate function into the model.

Thank you PerennialII. It's so nice to meet you and this forum. I'll come more often and update my progress.


Daniel

 

[PerennialII]

And forgot to welcome you to PF, belated one ! I'd say Comsol is build a bit different from other FE packages (pretty much all of them), and it doesn't "increment" or formulate the eq-systems like other software usually do (and I think that's one of the underlying issues affecting this as well, most software wouldn't complain about this, but then again they would lack many of the properties of Comsol). I've usually solved it simply by computing the derivatives like m-file(x+delta)-m-file(x-delta)/2*delta for models which are (for example) nonlinear and can't really get an analytical one (setting delta appropriately). And have gotten it to work with some pretty "weird" models.

Then again there are other tricks, developing things to series, adding stabilizing functions to aid differentiability etc., but those aren't really at the core of this problem.

 

[TWNien]

Hello, it's me again.
First, thanks for your hints. They are very helpful. So far I've made a big progress, and I want to make some remarks for this case. Actually, only one comment: the format of input and output. I was wrong that I thought COMSOL sends concentration of each points while calling m-function. In fact COMSOL sends c in the form of vector [c1 c2 c3 ... cn], that is I have to rewrite my m-function to deal with that.
After all things done, when I look back, it's kind of stupid stuck there for so long...

The next step will be using the better way to differentiate the m-function instead of using small delta.

 

[PerennialII]

Hello, it's me again.
First, thanks for your hints. They are very helpful. So far I've made a big progress, and I want to make some remarks for this case. Actually, only one comment: the format of input and output. I was wrong that I thought COMSOL sends concentration of each points while calling m-function. In fact COMSOL sends c in the form of vector [c1 c2 c3 ... cn], that is I have to rewrite my m-function to deal with that.
After all things done, when I look back, it's kind of stupid stuck there for so long...

The next step will be using the better way to differentiate the m-function instead of using small delta.

Yeah, the m-files should be made able to handle any form of input. Often get not-so-obvious error messages when fail to follow this "rule". Was doing some analyses of a highly nonlinear material last week and after a while of trying got it working best with a simple central - difference approach in an m-file combined with various smoothing functions to get it to accept a whole range of input. Although still favor developing pretty much everything complex in the m-files to a series (locally in a point of evaluation) and then working with them .