Integration over user defined function in GAUSS

The simplest method would be to use function intquad1 to calculate a one dimensional integral for example:

range_end = 1;
range_start = 0;
area = intquad1(&myProc, range_end|range_start);

proc (1) = myProc(x);
   retp(x .* sin(x));
endp;

The other option is to use inthp4. Here is an example. You can keep the wrapper function simpleInthp41d if you want and just call that function to make it simpler to use.

area = simpleInthp1D(&myProc, 0, 1);

proc (1) = myProc(struct DS *data, x);
   //If your function needs data other than from 'x'
      //you can pull it from the 'data' struct like this
      //new_var = data->dataMatrix;
   retp(x .* sin(x));
endp;

proc (1) = simpleInthp1D(fct, range_start, range_end);
   local fct:proc, area_;
	
   //Declare control structure and fill
   //in defaults
   struct inthpControl ctl;
   ctl = inthpControlCreate();

   //Make any control structure changes here
   //For example:
   //ctl.eps = 1e-5;
	
   //Declare and DS structure for any extra data
   struct DS d0;
   d0 = dsCreate();
	
   lim = range_end|range_start;
	
   //Call integration routine
   area_ = inthp4(&fct, &d0, ctl, lim);
   retp(area_);
endp;

Hello Please correct me if I am wrong, but in your examples isn't the integration is being done over a rectangular region. For example in your first example isn't it the assumption that x is uniformly distributed between 0 and 1?

What I want to do is I want to define a distribution for x other than rectangular, there are are options to do that for two/three dimensional integrals using intquad2 and intquad3, but I do not know how to use that in one dimenssinal case.

Look forward to your help.

Regards

Annesha

I am confused. When I first read your question, I was assuming that you simply wanted to calculate the area under the curve of a user defined function. For example, let's say the area under the curve of the sin function between 0 and 1.

Are you wanting to perform a double integral over a non-rectangular region as illustrated in slide number 2 at this link? http://www.math.ust.hk/~malwu/math100/lecture15.pdf

If not, could you point me to something that explains what you are trying to do?

Hello,

I do not want to do a double integral.

Let me try to explain what I want to do:

If g(.) is a function of x and f(.) is the pdf (probability density function of x), then all I want to do is to calculate the expected value of the function g(.) with respect to x, and the formula for doing this is

And also as I mentioned previously, I want to define both g(x) and the pdf, f(x). Hope its a little less ambiguous now.

I would appreciate your help in this regard.

Regards

Annesha

Finding expected value of a continuous function of a random variable
ExpectedValueRVFunc that works through an example analytically using the a uniform distribution function explicitly. Below is the code that corresponds to the analytical solution from the attachment.

//Main code
area = inthp4Simple(&myProc, 0, 10);
print area;

//User define procedures
proc (1) = pdfu(x, range_start, range_end);
   retp(1./(range_end - range_start));
endp;

proc (1) = myProc(struct DS *data, x);
   local range_end, range_start;
   //If your function needs data other than from 'x'
   //you can pull it from the 'data' struct like this
   range_start = data->dataMatrix[1];
   range_end = data->dataMatrix[2];
   retp(x .* x .* pdfu(x, range_start, range_end));
endp;

//Main helper procedure for reuse
proc (1) = inthp4Simple(fct, range_start, range_end);
   local fct:proc, area_, lim;
	
   //Declare control structure and fill
   //in defaults
   struct inthpControl ctl;
   ctl = inthpControlCreate();

   //Make any control structure changes here
   //For example:
   //ctl.eps = 1e-5;

   //Declare and DS structure which will
   //pass extra data to 'pdfu' function
   struct DS d0;
   d0 = dsCreate();
   //Passing range so it can be sent to 'pdfu'
   d0.dataMatrix = range_start|range_end;
	
   lim = range_end|range_start;

   //Call integration routine
   area_ = inthp4(&fct, &d0, ctl, lim);
   retp(area_);
endp;

If we make a couple of changes to the code above, we can make it work for an infinite integral and a non-uniform distribution. Here is an example using the normal distribution. You would replace pdfn with your probability density function:

//Main code
area = inthp4Simple(&myProc);
print area;

//User define procedure
proc (1) = myProc(struct DS *data, x);
   retp(.*pdfn(x));
endp;

//Main helper procedure for reuse
proc (1) = inthp4Simple(fct);
   local fct:proc, area_, lim;
	
   //Declare control structure and fill
   //in defaults
   struct inthpControl ctl;
   ctl = inthpControlCreate();
	
   //Declare and DS structure which will
   //pass extra data to 'pdfu' function
   struct DS d0;
   d0 = dsCreate();

   //Call integration routine
   //'inthp1' calculates an integral over
   //an infinite interval
   area_ = inthp1(&fct, &d0, ctl);
retp(area_);
endp;