1 function onvec = model_orthonormalize_gram_schmidt(model, model_data, vec,epsilon)
2 %
function onvec = model_orthonormalize_gram_schmidt(model, model_data, vec[,epsilon])
4 % Gram-Schmidt orthonormalization of vectors in columns of
5 % matrix vec to onvec. Almost zero vectors are deleted.
7 % note: a threshold epsilon is involved in
this routine
for detecting
8 % zero columns. This is a quite sensitive quantity! Change with care, as
9 % many functions build on
this!! Default is 1e-7
11 % Bernard Haasdonk 13.6.2002
12 % Martin Drohmann 14.5.2009
13 warning(
'gram-schmidt orthonormalization might be inaccurate!');
16 epsilon = 1e-7; % =>
for nonlinear evolution
this value is required.
18 %epsilon = 1e-10; =>
for linear evolution
this value was OK.
19 %epsilon = 1e-12; -> This causes error by returning non-orthogonal
20 % vectors (e.g. two identical vectors)
21 % So take epsilon larger than 1e-12!
24 onvec = zeros(size(vec));
28 A = model.get_inner_product_matrix(model_data);
30 % check on identity of vectors (can happen, that numerics does not detect
31 % this afterwards due to rounding errors!!)
33 for i=1:(size(vec,2)-1)
34 for j=(i+1):size(vec,2)
35 if isequal(vec(:,i),vec(:,j))
43 for i = 1:size(vec,2);
44 % orthogonalize next vector i and assume, that it is already
45 % orthogonal to previous ones
46 n = sqrt(onvec(:,i)'*A* onvec(:,i));
50 onvec(:,i) = onvec(:,i)/n;
53 % orthogonalize remaining vectors wrt this one:
55 A_mult_onvec_i = model.get_inner_product_matrix(model_data)*onvec(:,i);
57 % create row-vector with projections on the created on-vector
58 coeffs= A_mult_onvec_i' * onvec(:,i+1:end);
60 % create matrix of scaled versions of actual orthonormalized vector
61 coeffmat = onvec(:,i) * coeffs;
63 % perform orthogonalization
64 onvec(:,i+1:end) = onvec(:,i+1:end) - coeffmat;
68 % eliminate zero-columns
73 %resort zero-columns to back
74 %nsqr = sum(onvec.^2);
79 if max(max(abs(eye(size(onvec,2))-K)))>1e-3
81 error(['check orthonormalization!! Non orthogonal vectors' ...