Discussion Forum

Hi, man

Thanks a lot for sharing this well explained Matlab code!

One question: is it possible to do geometry parameter sweep with the code?

Thanks a lot for sharing this well explained Matlab code!

No problem, the question arises almost daily on the forum.

One question: is it possible to do geometry parameter sweep with the code?

Why not? You just have to redefine the geometry parameters and redo the meshing on every loop after the mphload command.

Hi Francois,

I'm a beginner in matalb with comsol.

I'd like to use a probe result as an entry for one of my material property. I used your matlab program and I adapted it. I run the model once, extract my variable 'dom1' (domain probe) in a vector E as
E(1)=mpheval(model,'dom1')
and then, I run the loop.

Unfortunately, Matlab returns an error 'no method 'set' with match signature found for class 'com.comsol.model.impl.MaterialModelImpl'.
E appears to be a structure and not a vector. I have no idea how I could use correctly the probe results...

I guess either it's me who did something wrong or i won't be able to do it this way...

I'll be glad if you could have a guess on that, or anyone else!

Hi Carole,

The 'mpheval' command returns a structure, not a vector. The vector, or the values that you are looking for, are within the structure created with 'mpheval'. With 'E(1)=mpheval(model,'dom1')' the structure 'E' will be created, within which there is the vector 'd1', the values that you want to assign to your material are within this vector, you retrieve them with E.d1(...).

But on your first loop, there is nothing in 'E' since it's a calculation from the results. You have to make sure that you initialize 'E' for the first loop...

Run a simple code (my example at the beginning of this thread) and see the structures created by 'mpheval' and play with it in order to see what's inside and how to retrieve stuff from it.

Cheers

Indeed, can't have carrots matching potatoes... :)
and Comsol's working with strings, while I was thinking numbers.

Thanks for the fast answer! my first comsol-loop on matlab is running!

Cheers

Hi Francois,

Thanks for the example. really helpful.

Still, i have problem with the comsol-matlab livelink. I used solid mechanics physic. I want to create a loop between comsol and matlab for every time step, that is by making every run equal to one time step. So, if I want to run the study for 1 second with 0.1 time step, then I will loop the program for 10 times.

My problem is, I can't update the value for the next time step since I don't know what code to be used.
Really appreciate if you or anyone know and willing to share on how to do it =)

Hi Khairul,

I guess you run stationary simulations? What values do you want to update on each run, boundary conditions? If yes, just look in the .m code of your model and you will easily find the lines where the values are assigned. Just retrieve data with mpheval at the end of a loop and assign new boundary conditions for the next loop.

Good luck

Hi Francois,
thanks for the reply.

I try to run time-dependent study actually. And the value that i want to update during each run is the initial value of its physics, that is for my case, the 'u' and 'du/dt'. I already did as per your advice and somehow help me a bit. But still I don't know/find the right matlab code to update the initial 'u' and 'du/dt'.

Again, thanks for your helpful advice.

Here are my thoughts:

1st option: If you have a 1D model, it's easy to polyfit 'u' and 'du/dt' as a function of x and then apply those functions for initial conditions. I have never done 2 or 3D polyfit, but I guess it's feasible with Matlab, you will definitely find the answer for that on the Matlab forum.

2nd option: Specifically applying values on each node of your domain. With 'mpheval', we retrieve every node values, but is it possible to manually apply those values on every nodes as initial condition for the next simulation? My guess is yes, but I never did this. In the Comsol GUI, I'm pretty sure that it's possible to set up the solver to use previous simulation results as initial conditions for a new run. Maybe you can set it up in the GUI and you won't have much to program in Matlab except that it's looping as long as you want. I hope someone else could elaborate on that, or there are probably good hints in other threads.

Yours,
François

Thanks a lot for sharing this well explained Matlab code!

No problem, the question arises almost daily on the forum.

One question: is it possible to do geometry parameter sweep with the code?

Why not? You just have to redefine the geometry parameters and redo the meshing on every loop after the mphload command.

Can you show me how to do geometry parameter sweep in MATLAB and redo the meshing?
When I refine the geometry parameters, use model.mesh('mesh1').run; it still cannot run and shows errors.

Hello,

Thanks for your answer to all the questions that have been asked from you. I don't have the livelink feature of COMSOL. So I am more trying to learn here. I want to know whether this feature help me in my problem.
Here is my problem:
I want to use COMSOL for one of my models. In this model, I want to find the viscosity and velocity in a geometry which I have both soil and water. So I should use Darcy's law. In this law velocity is dependent to viscosity. However, the viscosity which I need to find at the end is also dependent to Darcy's law velocity. Hence my fluid is non-Newtonian fluid.
When I add Darcy's law to my physics. It asked for fluid and matrix properties. For matrix properties I have no problems, but for fluid properties it asked for viscosity, which I should write in that spot a function of Darcy's law velocity ( f(darcy's law velocity) ). When I run this model the error is : ''circular variable dependency detected'. So I used the iterative solver instead of direct one. Still the same error. Does this livelink feature help me here? What I need to know is that in your code, you just enter ' k = {'200[W/m/K]' '1000[W/m/K]' '2000[W/m/K]'};' However my viscosity here is a function of Darcy's law velocity. How can I code this in MATLAB.
FYI, I am not professional in using MATLAB:(.

Thanks,

Negin

Hi Negin,

First, I think your question doesn't really fit in this thread. You should try to find a solution to your problem in the Comsol GUI first, then code it in Matlab afterwards.

Regarding your circular dependency problem, I think you should add a mathematic node (e.g. coefficient form PDE) that calculates the viscosity function of Darcy's velocity, and then in the Darcy's node you just refer to this variable to set the viscosity.

Thank you so much for your response. I will try your suggestion.

Negin

Hello, for completeness, here is what i found saving as M-file and later on accessing the COMSOL object directly.
I took the large deformation beam from the model library (structural mechanics verification models) and introduced a global parameter for the beam width, to loop over later in MATLAB, using the COMSOL GUI. The saved model M-file i modified to

function out = beamModel(W)
...
model.param.set('F_Lx', '-3.844[MN]', 'Maximum compressive load');
model.param.set('F_Ly', '1e-3*F_Lx', 'Transverse load');
model.param.set('NCL', '0', 'Normalized compressive load');
model.param.set('W', [num2str(W),'[N]'], 'beam width');
...

Then i run this model in a loop by another script (M-file):

clc; clear;
W=3.2; % this could be set in a loop
model=beamModel(W); % run beam model
% post-process using COMSOL commands (cf. LiveLink Docu)
mphplot(model,'pg1');
figure;
mphplot(model,'pg2');
ui=mphinterp(model,'u','coord',[W;0.005]); % 0.005 is half height
vi=mphinterp(model,'v','coord',[W;0.005]);
% alternatively mpheval for evaluation at the nodes
% post-process accessing the COMSOL object
model.result.numerical('pev1').set('expr','u');
ue=model.result.numerical('pev1').getReal(); % end point displacement (x) for each load step
model.result.numerical('pev1').set('expr','v');
ve=model.result.numerical('pev1').getReal(); % end point displacement (y) for each load step
figure; plot(ue, ve, ui, vi, 'Linewidth',2); % plot path of end point trajectory

I documented the steps more detailed for myself as tutorial. Maybe it is of interest for somebody else.

www.rz.uni-karlsruhe.de/~hf118/tutorial_comsol_matlab_livelink.pdf

Best Regards, D.

Hi Francois,


I have a simple model (coupled ode only). The model works fine in comsol and even in comsol-matlab interface.
The problem is when I try to extract my variable u using mpheval it gives error

??? Java exception occurred:
Exception:
java.lang.NullPointerException
(rethrown as com.comsol.util.exceptions.FlException)
Messages:



Stack trace:

at com.comsol.post.numerical.NumericalBuilder$PostEvalInterp.postCreateInit(Unknown
Source)

Can you please help me on this.

Thanks

parikshit

Hi Parikshit,

Look at the help for mpheval --> it's used for evaluating expressions on node points, and there is no nodes in your model!

You need to use mphglobal to retrieve your global quantities u and w.

Put the first and last lines in comments:
% function out = model
% out = model;

Then, after you ran your model, try this:
u = mphglobal(model,'mod1.u')

Thanks Francois! for a quick reply. It worked.

and I like the clarity with which you explained it.

thanks again!!!

Hi Parikshit

You shouldn't get a NullPointerException. I'll see to that Comsol returns a more descriptive error message in the future.


--
Lars Gregersen
Comsol Denmark

Hey, man.
Thanks for sharing this code.
It is durable to loop the material parameter like this. But this question bother me a long time, can you help me?
How to access into the time step loop and newton raphson loop from matlab code?

Let's say:
for t=1:t_total
for i=1:iter
...........
end iteration
end timestep


Thanks..

Xu

Hello -

Can you tell me where to put the new Matlab code/loop?
"k = {'200[W/m/K]' '1000[W/m/K]' '2000[W/m/K]'}; % strings to pass to Comsol
for n = 1:length(k)
model = mphload('example.mph'); % load model
ModelUtil.showProgress(true) % display the progress bar
model.material('mat1').propertyGroup('def').set('thermalconductivity', {k{n} '0' '0' '0' k{n} '0' '0' '0' k{n}}); % this is the line I noted in the 'preprocessing', it's modified on every loop
model.sol('sol1').runAll; % solve
T(n) = mpheval(model,'T'); % extract and store data
kxx(n) = mpheval(model,'ht.kxx'); % extract and store data
ModelUtil.remove('model') % remove model because it will be reloaded on every loop
end"

Should it be put at the bottom of example.m or within example.m or in another file?

Thank you!

Hi Rose,

You put it in the Command Window that opens when you start 'Comsol 4.2a with Matlab'. Make sure that 'example.mph' is in the current Matlab folder.

Francois -

Thank you! Is it possible to run it from a script?

Attached is the .mph corresponding .m and the loop I am trying to run (essentially save and analyze data through Matlab). There are lots of things in the run_code.m that could be done better, but any suggestions on taking the exported data, counting occurrences and then plotting would be useful (in this instance we are looking at the percentage of particles in each flow channel for different particle densities and diameters).

Well wishes and thank you!
Rosie

Thank you! Is it possible to run it from a script?

Definitely (and you should have tried before asking ;-)).

Attached is the .mph corresponding .m and the loop I am trying to run (essentially save and analyze data through Matlab). Any suggestions on what could be done better?

I don't have time for that, sorry.

Hi
I want to simulate bone remodeling in comsol. the question is how can I access time step and stress and strain from matlab code?

I would highly appreciate your help
Thanks, Mehran

Hi Francois,

Thanks a lot for this excellent post. It has been really useful for me.

Bruno Rango.

Hi Francois,
thanks for the reply.

I try to run time-dependent study actually. And the value that i want to update during each run is the initial value of its physics, that is for my case, the 'u' and 'du/dt'. I already did as per your advice and somehow help me a bit. But still I don't know/find the right matlab code to update the initial 'u' and 'du/dt'.

Again, thanks for your helpful advice.

Hello,

This posting is specifically to address the earlier question:
"I want to simulate bone remodeling in comsol. the question is how can I access time step and stress and strain from matlab code?"
With COMSOL version 5.2 there is now also the ability to define stress-strain relationships via external material functions coded in the C programming language. You may also want to consider this as an alternative approach.

Some resources to get you started:
www.comsol.com/blogs/accessing-external-material-models-for-structural-mechanics/
www.comsol.com/model/external-material-examples-structural-mechanics-32331

Best Regards,

Hi,
I am wondering when you reload the model, if the profile (velocity, temperature,pressure of all space) you calculated before would be loaded into the new loop? I was doing a 3D pipeflow and solid ht model. I need to use the outlet temperature (but not velocity) as the next inlet. So if i use a loop code like above, all the velocity, temperature will be back to what I set for the first loop. Is that right?
Thanks a lot.
Yang

Hi Francois
Thanks a lot for the post. It was ver helpful. I am actually working on the "Voltage induced in the coil by moving magnet" example given by comsol. In the default example, the magnet oscillates to produce voltage but for my project, I want to move the magnet inside the coil in a custom manner. I have the data of incremental displacements of the magnet and I stored the data in an array as u suggested in your post. I run the following code in the MATLAB workspace to change the displacements of the magnet-

k = {'5[cm]' '10[cm]' '12[cm]' '15[cm]' '12[cm]' '10[cm]' '5[cm]'}; % strings to pass to Comsol
for n = 1:length(k)
model = mphload('induced_voltage_moving_magnet.mph'); % load model
ModelUtil.showProgress(true) % display the progress bar
model.physics('ale').feature('disp2').setIndex('dx', 'k(n)', 1); % this is the line I noted in the 'preprocessing', it's modified on every loop
model.sol('sol1').runAll; % solve
ModelUtil.remove('model') % remove model because it will be reloaded on every loop
end

but I get the following result and error which I cannot understand. Plz help me resolve it.
ans =

1

Java exception occurred:
Exception:
com.comsol.util.exceptions.FlException: The following feature has encountered a
problem
(rethrown as com.comsol.util.exceptions.FlException)
(rethrown as com.comsol.util.exceptions.FlException)
Messages:
The following feature has encountered a problem

The following feature has encountered a problem

The following feature has encountered a problem:

Unknown function or operator.
- Name: k
- Feature: Stationary Solver 1 (sol1/s1)


Stack trace:

at com.comsol.clientapi.engine.c.handleException(Unknown Source)

at com.comsol.client.interfaces.f$d.f(Unknown Source)

at com.comsol.client.interfaces.f.a(Unknown Source)

at com.comsol.client.interfaces.f.runAndWait(Unknown Source)

at com.comsol.clientapi.engine.APIEngine.runMethod(Unknown Source)

at com.comsol.clientapi.engine.APIEngine.runMethod(Unknown Source)

at com.comsol.clientapi.impl.SolverSequenceClient.runAll(Unknown Source)

Caused by: Exception:
com.comsol.util.exceptions.FlException: The following feature has encountered a
problem
(rethrown as com.comsol.util.exceptions.FlException)
Messages:
The following feature has encountered a problem

The following feature has encountered a problem:

Unknown function or operator.
- Name: k
- Feature: Stationary Solver 1 (sol1/s1)


at com.comsol.clientapi.engine.c.handleException(Unknown Source)

at com.comsol.client.interfaces.f.a(Unknown Source)

at com.comsol.client.interfaces.f.processCommandAnswer(Unknown Source)

... 6 more

Caused by: Exception:
com.comsol.util.exceptions.FlException: The following feature has encountered a
problem
Messages:
The following feature has encountered a problem:

Unknown function or operator.
- Name: k
- Feature: Stationary Solver 1 (sol1/s1)


at com.comsol.solver.SolverOperation.addError(Unknown Source)

at com.comsol.solver.SolverOperation.execute(Unknown Source)

at com.comsol.model.method.SolverSequenceMethod.doRun(Unknown Source)

at com.comsol.model.internal.impl.SolverSequenceImpl.a(Unknown Source)

at com.comsol.model.internal.impl.SolverSequenceImpl.k(Unknown Source)

at com.comsol.model.internal.impl.SolverSequenceImpl$23.a(Unknown Source)

at com.comsol.model.internal.impl.SolverSequenceImpl$23.execute(Unknown Source)

at com.comsol.model.clientserver.ClientManager$1.call(Unknown Source)

at java.util.concurrent.FutureTask.run(Unknown Source)

at java.util.concurrent.ThreadPoolExecutor.runWorker(Unknown Source)

at java.util.concurrent.ThreadPoolExecutor$Worker.run(Unknown Source)

at java.lang.Thread.run(Unknown Source)

Thanks a lot for the post.

You are welcome!

model.physics('ale').feature('disp2').setIndex('dx', 'k(n)', 1); % this is the line I noted in the 'preprocessing',

Maybe the error is because you wrote k(n) instead of k{n} ?

I don't use Comsol anymore, never used version 5.x, my Comsol 'experience' ended at version 4.2a.

François