Most of these look sensible. I'd be inclined to try VEGAS and compare against plain Monte Carlo.
On Mon, Jan 1, 2018 at 6:10 AM, Vasu Jaganath <[email protected]> wrote: > HI Martin, Yes one of my Q is very discontinuous with respect to my > integration variable Z. > > I have attached the plots for you to see > > On Sun, Dec 31, 2017 at 9:20 PM, Vasu Jaganath <[email protected]> > wrote: > >> Yes, I will show you the plots soon, Q is actually 2 variable function >> but for my calculations I am treating one of the variables as a parameter, >> which is a physically valid assumption. Yes I do encounter alpha and beta >> values less than 1. >> >> On Sun, Dec 31, 2017, 9:13 PM Martin Jansche <[email protected]> wrote: >> >>> So you want to find E[f] = \int_0^1 f(x) dbeta(x | a, b) dx. Can you plot >>> your typical f? And what are typical values of the parameters a and b? Do >>> you encounter a<=1 or b<=1? If so, how does f(x) behave as x approaches 0 >>> or 1? >>> >>> On Mon, Jan 1, 2018 at 3:37 AM, Patrick Alken <[email protected]> >>> wrote: >>> >>> > The question is whether your Q contains any singularities, or is highly >>> > oscillatory? Is such cases fixed point quadrature generally doesn't do >>> > well. If Q varies fairly smoothly over your interval, you should give >>> > fixed point quadrature a try and report back if it works well enough >>> for >>> > your problem. The routines you want are documented here: >>> > >>> > http://www.gnu.org/software/gsl/doc/html/integration.html# >>> > fixed-point-quadratures >>> > >>> > Also, if QAGS isn't working well for you, try also the CQUAD routines. >>> > I've found CQUAD is more robust than QAGS in some cases >>> > >>> > On 12/31/2017 05:28 PM, Vasu Jaganath wrote: >>> > > I have attached my entire betaIntegrand function. It is a bit >>> complicated >>> > > and very boiler-plate, It's OpenFOAM code (where scalar = double >>> etc.) I >>> > > hope you can get the jist from it. >>> > > I can integrate the Q using monte-carlo sampling integration. >>> > > >>> > > Q is nothing but tabulated values of p,rho, mu etc. I lookup Q using >>> the >>> > > object "solver" in the snippet. >>> > > >>> > > I have verified evaluating <Q> when I am not using those <Q> values >>> back >>> > in >>> > > the solution, It works OK. >>> > > >>> > > Please ask me anything if it seems unclear. >>> > > >>> > > >>> > > >>> > > >>> > > >>> > > >>> > > On Sun, Dec 31, 2017 at 3:55 PM, Martin Jansche <[email protected]> >>> > wrote: >>> > > >>> > >> Can you give a concrete example of a typical function Q? >>> > >> >>> > >> On Sat, Dec 30, 2017 at 3:42 AM, Vasu Jaganath < >>> > [email protected]> >>> > >> wrote: >>> > >> >>> > >>> Hi forum, >>> > >>> >>> > >>> I am trying to integrate moments, basically first order moments >>> <Q>, >>> > i.e. >>> > >>> averages of some flow fields like temperature, density and mu. I am >>> > >>> assuming they distributed according to beta-PDF which is >>> parameterized >>> > on >>> > >>> variable Z, whose mean and variance i am calculating separately and >>> > using >>> > >>> it to define the shape of the beta-PDF, I have a cut off of not >>> using >>> > the >>> > >>> beta-PDF when my mean Z value, i.e <Z> is less than a threshold. >>> > >>> >>> > >>> I am using qags, the adaptive integration routine to calculate the >>> > moment >>> > >>> integral, however I am restricted to threshold of <Z> = 1e-2. >>> > >>> >>> > >>> It complains that the integral is too slowly convergent. However >>> > >>> physically >>> > >>> my threshold should be around 5e-5 atleast. >>> > >>> >>> > >>> I can integrate these moments with threshold upto 5e-6, using >>> > Monte-Carlo >>> > >>> integration, by generating random numbers which are >>> beta-distributed. >>> > >>> >>> > >>> Should I look into fixed point integration routines? What routines >>> > would >>> > >>> you suggest? >>> > >>> >>> > >>> Here is the (very simplified) code snippet where, I calculate >>> alpha and >>> > >>> beta parameter of the PDF >>> > >>> >>> > >>> zvar = min(zvar,0.9999*zvar_lim); >>> > >>> alpha = zmean*((zmean*(1 - zmean)/zvar) - 1); >>> > >>> beta = (1 - zmean)*alpha/zmean; >>> > >>> >>> > >>> // inside the fucntion to be integrated >>> > >>> // lots of boilerplate for Q(x) >>> > >>> f = Q(x) * gsl_ran_beta_pdf(x, alpha, beta); >>> > >>> >>> > >>> // my integration call >>> > >>> >>> > >>> helper::gsl_integration_qags (&F, 0, 1, 0, 1e-2, >>> > 1000, >>> > >>> w, &result, >>> &error); >>> > >>> >>> > >>> And also, I had to give relative error pretty large, 1e-2. However >>> > some of >>> > >>> Qs are pretty big order of 1e6. >>> > >>> >>> > >>> Thanks, >>> > >>> Vasu >>> > >>> >>> > >> >>> > >>> > >>> > >>> >> >
