Showing results for

- Intel Community
- Software Development SDKs and Libraries
- Intel® oneAPI Math Kernel Library & Intel® Math Kernel Library
- GaussianMV parameter error

- Subscribe to RSS Feed
- Mark Topic as New
- Mark Topic as Read
- Float this Topic for Current User
- Bookmark
- Subscribe
- Mute
- Printer Friendly Page

Highlighted
##

Hello, I am trying to work with the random number generator of the mutlivariate Gaussian distribution in FORTAN. TO this end, I have compiled my MKL library, and am attempting to get the random numbers.

1. I have several questions. One it is hard to tell from the documentation. whether ytou need the Cholesky decompostion, or the full covariance matrix. which is required

2. When I attempt to run the below program, I recive the following error: Parameter 4 was incorrect on entry to vsRngGaussia

Does anyone have any suggestions as to what is wrong with my code. Please also note that I have tried using the Cholesky decomposition code from LAPACK but it does not seem to work.

jaureus

Beginner

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Email to a Friend
- Report Inappropriate Content

04-25-2011
02:19 PM

7 Views

GaussianMV parameter error

1. I have several questions. One it is hard to tell from the documentation. whether ytou need the Cholesky decompostion, or the full covariance matrix. which is required

2. When I attempt to run the below program, I recive the following error: Parameter 4 was incorrect on entry to vsRngGaussia

Does anyone have any suggestions as to what is wrong with my code. Please also note that I have tried using the Cholesky decomposition code from LAPACK but it does not seem to work.

program mcnet

implicit none

integer :: i, j

real :: harvest1, harvest2, sd

integer :: dimen, seed_in1, seed_in2, n_rand_obs_in

real , dimension (:, :), allocatable :: newidentity, sigmanew, meanvectnorm, meanvectnormtrans, meanvect, meanvecttrans, &

meanvectold, meanvecttransold, Q, r_mvgaussian_list1, r_mvgaussian_list2, sigmanewchol, newidentity2

interface

subroutine gaussianmv_rng_mkl(seed_in, n_rand_obs_in, r_mvgaussian_list, dimen, mean, tcv)

implicit none

integer, intent(in) :: seed_in, n_rand_obs_in, dimen

real, dimension(1:n_rand_obs_in, 1:dimen), intent(out) :: r_mvgaussian_list

real, dimension(1:1, 1:dimen), intent(in) :: mean

real, dimension(1:dimen, 1:dimen), intent(in) :: tcv

end subroutine gaussianmv_rng_mkl

end interface

! allowing for easy changes in dimensionality

dimen = 3

n_rand_obs_in = 3

allocate(newidentity(1:dimen, 1:dimen))

!allocate(sigmanew(1:dimen, 1:dimen))

allocate(meanvectnorm(1:dimen, 1:1))

allocate(meanvectnormtrans(1:1, 1:dimen))

!allocate(meanvect(1:dimen, 1:1))

!allocate(meanvecttrans(1:1, 1:dimen))

!allocate(meanvectold(1:dimen, 1:1))

!allocate(meanvecttransold(1:1, 1:dimen))

!allocate(sigmanewchol(1:dimen, 1:dimen))

allocate(newidentity2(1:dimen, 1:dimen))

allocate(r_mvgaussian_list2(1:n_rand_obs_in, 1:dimen) )

allocate(Q(1:dimen, 1:1))

!allocate(tcv(1:dimen, 1:dimen))

newidentity = 1.0

!scaling factor as per algoritm

sd = (2.38)**2/dimen

meanvect = 0.0

meanvectnorm = 0.0

do j = 1, dimen

! identity matrix times scaling factor as per Roberts and Rosthenthal (2009) gives innitial value of covariance matrix

newidentity( j, j) = 1.0*(0.01)/dimen

end do

print *, newidentity

call random_number(harvest1)

seed_in1 = aint(harvest1 * 2.0e0**31)

meanvectnorm(1,1) = 0.2

meanvectnorm(2,1) = 0.3

meanvectnorm(3,1) = 0.001

!print *, meanvectnorm

meanvectnormtrans = transpose(meanvectnorm)

! call spotrf('U', dimen, newidentity)

!call potrf( newidentity, 'U')

! note that since cholesky decomposition of newidentity does not change, since newidentity does not change one only needs to derive it once

newidentity2 = 0.0

do i = 1, dimen

do j = 1, dimen

if (i .ge. j) then

newidentity2(i,j) = newidentity(i,j)

end if

end do

end do

print *, newidentity2

print *, "test"

call gaussianmv_rng_mkl(seed_in1, n_rand_obs_in, r_mvgaussian_list2, dimen, meanvectnormtrans, newidentity2)

!subroutine gaussianmv_rng_mkl(seed_in, n_rand_obs_in, r_mvgaussian_list, dimen, mean, tcv)

do i = 1, n_rand_obs_in

do j = 1, dimen

print *, r_mvgaussian_list2(i, j)

end do

end do

end program

which in turns calls the subroutine that I have created to use mutlivariate random numbers,

subroutine gaussianmv_rng_mkl(seed_in, n_rand_obs_in, r_mvgaussian_list, dimen, mean, tcv)

use :: mkl_vsl_type

integer :: brng

integer, intent(in) :: seed_in, n_rand_obs_in

real, pointer, intent(out) :: r_mvgaussian_list

type(VSL_STREAM_STATE) :: stream

integer :: status, methodrng

external :: vslnewstream, vsrnggausianmv, vsldeletestream

methodrng = VSL_RNG_METHOD_GAUSSIANMV_BOXMULLER

mstorage = VSL_MATRIX_STORAGE_FULL

brng = vsl_brng_mcg31

status = vslnewstream(stream, brng, seed_in)

status = vsrnggaussianmv( methodrng, stream, n_rand_obs_in, r_mvgaussian_list, dimen, mstorage, mean, tcv)

status = vsldeletestream(stream)

end subroutine

2 Replies

Highlighted
##

yuriisig

Beginner

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Email to a Friend
- Report Inappropriate Content

04-26-2011
01:01 AM

7 Views

see ....\mkl\examples\vslf

Highlighted
##

You have a large number of errors with respect to argument types in your sample code. To catch such errors at compile time rather than experiencing run-time aborts, use IMPLICIT NONE, use all available interface modules, and use /warn:all as a compiler option. Finally, note that one of the worst things you can do is to give the compiler an incorrect interface block.

The documentation clearly states that you need to supply the Cholesky factor T of the variance-covariance matrix as an argument.

There is also a functioning example,*vsrnggaussianmv.f*, provided in the MKL examples directory. The example shows the call the call to the Cholesky factorization routine *spotrf*.

mecej4

Black Belt

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Email to a Friend
- Report Inappropriate Content

04-26-2011
01:06 AM

7 Views

The documentation clearly states that you need to supply the Cholesky factor T of the variance-covariance matrix as an argument.

There is also a functioning example,

For more complete information about compiler optimizations, see our Optimization Notice.