The Origin Forum
File Exchange
Try Origin for Free
The Origin Forum
Home | Profile | Register | Active Topics | Members | Search | FAQ | Send File to Tech support
Username:
Password:
Save Password
Forgot your Password? | Admin Options

 All Forums
 Origin Forum
 Origin Forum
 Unexpected diverging integral using INTEGRAL
 New Topic  Reply to Topic
 Printer Friendly
Author Previous Topic Topic Next Topic Lock Topic Edit Topic Delete Topic New Topic Reply to Topic

johnnyho1415

Hong Kong
3 Posts

Posted - 06/19/2019 :  12:38:29 AM  Show Profile  Edit Topic  Reply with Quote  View user's IP address  Delete Topic
Origin Ver. and Service Release (OriginPro 2018 (32-bit) SR1, b9.5.1.195 / OriginPro 2017 (64-bit) SR2, b9.4.2.380):
Operating System: Win7 Pro (32-bit) / Win7 Home (64-bit)

Something funny happens if you try this LabTalk script:

wks.ncols = 4;
col(a)={1.5,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20,25,30,35,40,45,50,60,70,80,90,100,200,300,400,500,1000,2000,3000,4000,5000};
wks.col3.type = 4;
col(c)={0:0.02:1.5};

%Z = "double k = 8.6173303E-5;%(CRLF)";
%Z = "%Zdouble Tcell = 300;%(CRLF)";
%Z = "%Zfunction double Fc1(double energy, double Eg, double kT) {%(CRLF)";
%Z = "%Z	double Absorptance = (energy-Eg)>0 ? 1 : 0;%(CRLF)";
%Z = "%Z	return (Absorptance==0||energy/kT>287*ln(10))?0:Absorptance*energy*energy/(exp(energy/kT)-1);%(CRLF)";
%Z = "%Z}";
csetvalue f:="Integral(Fc1, 0, col(a), 1.5, k*Tcell)" c:=2 r:=1 s:=%Z;
csetvalue f:="Integral(Fc1, col(c), +inf, 1.5, k*Tcell)" c:=4 r:=1 s:=%Z;

plotxy iy:=[Book1]1!(1,2) o:=[Graph1]1! plot:=202;
layer.x.type = 2;
layer.y.type = 2;
layer.x.labelSubtype = 2;
layer.y.labelSubtype = 2;
layer.y.inc = 10;
layer.y.minorTicks = 0;
Rescale;

plotxy iy:=[Book1]1!(3,4) o:=[Graph2]1! plot:=202;
layer.y.type = 2;
layer.y.labelSubtype = 2;
//layer.y.inc = 10;
layer.y.minorTicks = 0;
Rescale;


The non-negative function involved effectively ensembles the blackbody radiation (w/o constants) at 300 K defined at energy > 1.5 eV with the counts integrated. The function contains no singularities (except 0) and converges to 0 at +inf.

In the above script, Col(A) & Col(B) simulate the range of integration from 0 to a varying large number (to mimic +inf) to check for convergence. Col(C) & Col(D) simulate the range of integration from a varying small number (to mimic 0) where the function value is 0 to +inf.

It should be expected that both of them should converge to the number of the range (0,+inf) of the value around 3.808E-27. However, both turn out not to converge at all. I have checked that the function value is consistent, but the integral value is funny.

May I ask which algorithm the command
Integral
uses for evaluating the integral? Am I missing out something, or this command appears tricky? Thanks very much.


Johnny

Edited by - johnnyho1415 on 06/19/2019 12:41:33 AM

cpyang

USA
1406 Posts

Posted - 06/19/2019 :  8:03:41 PM  Show Profile  Edit Reply  Reply with Quote  View user's IP address  Delete Reply
Instead of all these code, maybe you can provide the data and show the Integral function result what you expect?

Can you paste the data here or email?

CP
Go to Top of Page

johnnyho1415

Hong Kong
3 Posts

Posted - 06/20/2019 :  03:12:54 AM  Show Profile  Edit Reply  Reply with Quote  View user's IP address  Delete Reply
Hi CP,

Thanks for your response. As I encountered issues on attaching files here, I chose to generate the data with the LabTalk scripts so that everyone can reproduce the problem easily without ambiguity.

You will have the worksheet of values right away when you run the first two blocks of the script in my first post (Start a new project -> Press Alt+Shift+3 -> Paste the script -> Select all (Ctrl+A) -> Press Enter):

wks.ncols = 4;
col(a)={1.5,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20,25,30,35,40,45,50,60,70,80,90,100,200,300,400,500,1000,2000,3000,4000,5000};
wks.col3.type = 4;
col(c)={0:0.02:1.5};

%Z = "double k = 8.6173303E-5;%(CRLF)";
%Z = "%Zdouble Tcell = 300;%(CRLF)";
%Z = "%Zfunction double Fc1(double energy, double Eg, double kT) {%(CRLF)";
%Z = "%Z	double Absorptance = (energy-Eg)>0 ? 1 : 0;%(CRLF)";
%Z = "%Z	return (Absorptance==0||energy/kT>287*ln(10))?0:Absorptance*energy*energy/(exp(energy/kT)-1);%(CRLF)";
%Z = "%Z}";
csetvalue f:="Integral(Fc1, 0, col(a), 1.5, k*Tcell)" c:=2 r:=1 s:=%Z;
csetvalue f:="Integral(Fc1, col(c), +inf, 1.5, k*Tcell)" c:=4 r:=1 s:=%Z;


You can see the formulae in the worksheet as usual with the Set Column Values (Ctrl+Q) menu. I know that it is more readable. Hope this relieves you a bit.

The last 2 blocks of script just generates the plots of integral values as a function of its upper/lower limit for better readability.

---------------------------

Being new here, is there any advice on attaching images in the posts? (I tried, the images can be uploaded, but invoking the images (preview mode) did not have the content displayed. They are png.) I hoped to show the expected results but couldn't because of this technical issue. I believe solving this could help idea exchange.

Thanks again!


Johnny
Go to Top of Page

cpyang

USA
1406 Posts

Posted - 06/20/2019 :  2:46:46 PM  Show Profile  Edit Reply  Reply with Quote  View user's IP address  Delete Reply
We added an internal support item

https://originlab.jira.com/browse/SUP-3191

and we will investigate. Will let you know shortly.

CP
Go to Top of Page

johnnyho1415

Hong Kong
3 Posts

Posted - 06/21/2019 :  07:27:24 AM  Show Profile  Edit Reply  Reply with Quote  View user's IP address  Delete Reply
Thanks CP. Do I need to get through an authentication process to access to your URL provided? Or I don't need to keep up with that?
Go to Top of Page

cpyang

USA
1406 Posts

Posted - 06/21/2019 :  2:03:02 PM  Show Profile  Edit Reply  Reply with Quote  View user's IP address  Delete Reply
Hi Johnny,

that link is really for our internal use, and we will get back to you once our team has a chance to study the issue.

CP
Go to Top of Page

Sam Fang

293 Posts

Posted - 06/23/2019 :  11:38:21 PM  Show Profile  Edit Reply  Reply with Quote  View user's IP address  Delete Reply
Hi Johnny,

It was caused by your discontinuous integral function. There are two discontinuous points (x=1.5 and x=k*Tcell*287*ln(10)) in your integral function.

Origin's integral function uses adaptive quadrature method using the Gauss 10-point and Kronrod 21-point rules, and sample points for quadrature are not uniform. If there are discontinuous points or a sharp peak in the integral function, it may fail to produce sample points in the specified interval around these points.

For your case, you can divide the integral interval into several sections to fix the discontinuous points issue.
You can first define the second discontinuous point in Before Formula Scripts:
double xc=k*Tcell*287*ln(10);


For col(B), you can update the formula as:
Integral(Fc1, 0, 1.5, 1.5, k*Tcell)+ 
 ( col(A)<xc? Integral(Fc1, 1.5, col(a), 1.5, k*Tcell) : 
  ( Integral(Fc1, 1.5, xc, 1.5, k*Tcell)+
    Integral(Fc1, xc, col(A), 1.5, k*Tcell) ) )


For col(D), you can update the formula as:
Integral(Fc1, col(c), 1.5, 1.5, k*Tcell)+
 Integral(Fc1, 1.5, xc, 1.5, k*Tcell)+
 Integral(Fc1, xc, +inf, 1.5, k*Tcell)


In this way, you will find both of converged integral results are:
3.85112E-27.

Sam
OriginLab Technical Services

Edited by - Sam Fang on 06/23/2019 11:46:41 PM
Go to Top of Page
  Previous Topic Topic Next Topic Lock Topic Edit Topic Delete Topic New Topic Reply to Topic
 New Topic  Reply to Topic
 Printer Friendly
Jump To:
The Origin Forum © 2020 Originlab Corporation Go To Top Of Page
Snitz Forums 2000