|
[Sponsors] |
November 9, 2015, 15:50 |
Integration(math) issue
|
#1 | |
Member
Sandra
Join Date: Oct 2014
Posts: 58
Rep Power: 11 |
Hello Foamers!
This code might look familiar but I have added a few lines that require an integration that isn't going so well. Anyhoo, here's the code Code:
forAll(x3, celli) { for (count=0; count<iter; count++) { x3 = (x1 + x2)/2; jc = i*(Foam::exp(8.1*x))-(Foam::exp(-2.7*x)) f1 = (Foam::exp(8.1*x1))-(Foam::exp(-2.7*x1))-(jc/i); f3 = (Foam::exp(8.1*x3))-(Foam::exp(-2.7*x3))-(jc/i); if(neg(f1[celli] *f3[celli])) { x2 = x3; } else if((mag(x1[celli] - x2[celli]) < 1e-7)||(f3[celli]==0)||(jc[celli].integrate(0,3.6e9)==2e4)) { root =x3; break; } else { x1= x3; } } Info<<"Root is"<<endl; Info<<root<<endl; } Quote:
Sandra |
||
November 10, 2015, 03:20 |
|
#2 |
Senior Member
|
Hi,
1. Code:
f3[celli]==0 Code:
inline bool equal(const Scalar& s1, const Scalar& s2) { return mag(s1 - s2) <= ScalarVSMALL; } inline bool notEqual(const Scalar s1, const Scalar s2) { return mag(s1 - s2) > ScalarVSMALL; } And again this Code:
jc[celli].integrate(0,3.6e9)==2e4 |
|
November 10, 2015, 03:51 |
|
#3 | |
Member
Sandra
Join Date: Oct 2014
Posts: 58
Rep Power: 11 |
Quote:
Lastly, out of curiosity, why would Code:
jc[celli].integrate(0,3.6e9)==2e4 |
||
November 10, 2015, 04:10 |
|
#4 | |
Senior Member
|
If, as you said, jc is volScalarField, then jc[celli] can not be "...not a linear function per say...". It is a number and on a given time step it is in general even constant. I guess, you need to explain what you are trying to achieve with jc[celli].integrate(...). Right now what you are trying to do is like this:
Code:
(3.1415).integrate(0, 3e9) Quote:
https://en.wikipedia.org/wiki/Floati...uracy_problems https://duckduckgo.com/?q=floating+point+comparison |
||
November 10, 2015, 04:27 |
|
#5 |
Member
Sandra
Join Date: Oct 2014
Posts: 58
Rep Power: 11 |
sorry, my bad....
jc in paraview (at each time step) is a curve like a straight up witch hat and I want to iterate until the area under that curve equals a constant that I give it (in this case, 20,000). So, I agree, jc[celli] wouldn't work here but I wasn't thinking of that, unfortunately, when I tried it after "jc.integrate" gave me an error along the lines of "member doesn't have function integrate" |
|
November 10, 2015, 04:37 |
|
#6 |
Senior Member
|
Well, let us start from the beginning...
What is jc? In the first post you wrote, that it is volScalarField, i.e. something that varies in VOLUME, then you have written, that it is "...a curve like a straight up witch hat...", i.e. someting that varies just in ONE direction (surely there are 3D curves, yet they usually can be reduced to function of one variable). After you define jc, I think, it is easier to invent a method for calculation your "integral". |
|
November 10, 2015, 04:41 |
|
#7 |
Member
Sandra
Join Date: Oct 2014
Posts: 58
Rep Power: 11 |
it is defined as a volScalarField and calculated from the equation in the initial post BUT my case is 2D ie I have front and back patches empty so technically it's an areaScalarField if that were a thing in OpenFOAM and that's why I talked about a 2d curve.
Does that help? Sandra |
|
November 10, 2015, 04:58 |
|
#8 |
Senior Member
|
OK, your case is 2D, so is jc a surface? So you have
and would like to calculate where xa, xb, ya, yb are the corners of your simulation domain (so it would be volume under the curve and in simplest form it is a sum of jc[celli]*mesh.V()[celli]/<distance in empty direction>)? Or you would like to take values of jc along the line L (technically reducing surface to curve) and them integrate 1D function along L? |
|
November 10, 2015, 05:03 |
|
#9 |
Member
Sandra
Join Date: Oct 2014
Posts: 58
Rep Power: 11 |
Quote:
Many thanks, BTW |
|
November 10, 2015, 05:30 |
|
#10 |
Senior Member
|
There are several possibilities:
1. Your mesh is static regular orthogonal, your line is parallel to axis. You construct cell set (cells, where you line passes) yourself in the beginning of the simulation, starting from xa (using, for example, findCell method of mesh) and then iterating through the neighbors and selecting cells in required direction (or if mesh is not large, you can use findCell for every point of your line, as you are doing it just once). After in simplest case you iterate through the labels, take value of jc in the cell, take width of the cell and multiply. 2. You mesh is static irregular, your line is defined as a*x + b. You can use polyLineSet set class to construct cell set. 3. You mesh moves, is irregular, your line is defined as a*x + b. You reconstruct cell set using polyLineSet after every mesh topological change. |
|
|
|
Similar Threads | ||||
Thread | Thread Starter | Forum | Replies | Last Post |
[swak4Foam] swak4Foam 0.3.1 OF2.2.2 Installation Issue | scnsc | OpenFOAM Community Contributions | 3 | August 20, 2014 07:59 |
Strange issue while launching ANSYS workbench in CentOS 6.4 | Philip_C | ANSYS | 11 | August 29, 2013 06:44 |
CyclicAMI Issue In OpenFOAM 2.2.0 | prasant | OpenFOAM Running, Solving & CFD | 17 | March 16, 2013 02:00 |
Pressure boundary condition issue | Vijay | FLUENT | 0 | April 6, 2012 13:35 |
Meshing related issue in Flow EFD | appu | FloEFD, FloWorks & FloTHERM | 1 | May 22, 2011 08:27 |