mirror of
https://github.com/revarbat/BOSL2.git
synced 2025-01-16 13:50:23 +01:00
fix tests and bugs revealed by new tests
This commit is contained in:
parent
a5ae4879be
commit
d66c1929a5
@ -545,10 +545,12 @@ function gaussian_rands(N=1, mean=0, cov=1, seed=undef) =
|
|||||||
rdata = [for (i = count(dim*N,0,2)) sqrt(-2*ln(nums[i]))*cos(360*nums[i+1])]
|
rdata = [for (i = count(dim*N,0,2)) sqrt(-2*ln(nums[i]))*cos(360*nums[i+1])]
|
||||||
)
|
)
|
||||||
dim==1 ? add_scalar(sqrt(cov)*rdata,mean) :
|
dim==1 ? add_scalar(sqrt(cov)*rdata,mean) :
|
||||||
|
assert(is_matrix_symmetric(cov),"Supplied covariance matrix is not symmetric")
|
||||||
let(
|
let(
|
||||||
L = cholesky(cov)
|
L = cholesky(cov)
|
||||||
)
|
)
|
||||||
array_group(rdata,dim)*transpose(L);
|
assert(is_def(L), "Supplied covariance matrix is not positive definite")
|
||||||
|
move(mean,array_group(rdata,dim)*transpose(L));
|
||||||
|
|
||||||
|
|
||||||
// Function: spherical_random_points()
|
// Function: spherical_random_points()
|
||||||
@ -1092,14 +1094,13 @@ function _back_substitute(R, b, x=[]) =
|
|||||||
function cholesky(A) =
|
function cholesky(A) =
|
||||||
assert(is_matrix(A,square=true),"A must be a square matrix")
|
assert(is_matrix(A,square=true),"A must be a square matrix")
|
||||||
assert(is_matrix_symmetric(A),"Cholesky factorization requires a symmetric matrix")
|
assert(is_matrix_symmetric(A),"Cholesky factorization requires a symmetric matrix")
|
||||||
echo(A=A,len=len(A))
|
|
||||||
_cholesky(A,ident(len(A)), len(A));
|
_cholesky(A,ident(len(A)), len(A));
|
||||||
|
|
||||||
function _cholesky(A,L,n) = let(ffee=echo(insideA=A,L,n))
|
function _cholesky(A,L,n) =
|
||||||
A[0][0]<0 ? undef : // Matrix not positive definite
|
A[0][0]<0 ? undef : // Matrix not positive definite
|
||||||
len(A) == 1 ? submatrix_set(L,[[sqrt(A[0][0])]], n-1,n-1):
|
len(A) == 1 ? submatrix_set(L,[[sqrt(A[0][0])]], n-1,n-1):
|
||||||
let(
|
let(
|
||||||
i = n+1-len(A),ff=echo(i=i,lenA=len(A))
|
i = n+1-len(A)
|
||||||
)
|
)
|
||||||
let(
|
let(
|
||||||
sqrtAii = sqrt(A[0][0]),
|
sqrtAii = sqrt(A[0][0]),
|
||||||
|
@ -368,6 +368,17 @@ module test_gaussian_rands() {
|
|||||||
assert_equal(len(nums3), 1000);
|
assert_equal(len(nums3), 1000);
|
||||||
assert_equal(nums1, nums3);
|
assert_equal(nums1, nums3);
|
||||||
assert(nums1!=nums2);
|
assert(nums1!=nums2);
|
||||||
|
|
||||||
|
R = [[4,2],[2,17]];
|
||||||
|
data = gaussian_rands(100000,[0,0],R,seed=49);
|
||||||
|
assert(approx(mean(data), [0,0], eps=1e-2));
|
||||||
|
assert(approx(transpose(data)*data/len(data), R, eps=2e-2));
|
||||||
|
|
||||||
|
R2 = [[4,2,-1],[2,17,4],[-1,4,11]];
|
||||||
|
data3 = gaussian_rands(100000,[1,2,3],R2,seed=97);
|
||||||
|
assert(approx(mean(data3),[1,2,3], eps=1e-2));
|
||||||
|
cdata = move(-mean(data3),data3);
|
||||||
|
assert(approx(transpose(cdata)*cdata/len(cdata),R2,eps=.1));
|
||||||
}
|
}
|
||||||
test_gaussian_rands();
|
test_gaussian_rands();
|
||||||
|
|
||||||
|
Loading…
x
Reference in New Issue
Block a user