| 108 | |
| 109 | |
| 110 | UD_scheme::Float |
| 111 | UD_scheme::predictGq (const Matrix& Fx, const Matrix& G, const FM::Vec& q) |
| 112 | /* MWG-S prediction from Bierman p.132 |
| 113 | * q can have order less then x and a matching G so GqG' has order of x |
| 114 | * Precond: |
| 115 | * UD |
| 116 | * Postcond: |
| 117 | * UD |
| 118 | * |
| 119 | * Return: |
| 120 | * reciprocal condition number, -1 if negative, 0 if semi-definite (including zero) |
| 121 | */ |
| 122 | { |
| 123 | std::size_t i,j,k; |
| 124 | const std::size_t n = x.size(); |
| 125 | const std::size_t Nq = q.size(); |
| 126 | const std::size_t N = n+Nq; |
| 127 | Float e; |
| 128 | // Check preallocated space for q size |
| 129 | if (Nq > q_max) |
| 130 | error (Logic_exception("Predict model q larger than preallocated space")); |
| 131 | |
| 132 | if (n > 0) // Simplify reverse loop termination |
| 133 | { |
| 134 | // Augment d with q, UD with G |
| 135 | for (i = 0; i < Nq; ++i) // 0..Nq-1 |
| 136 | { |
| 137 | d[i+n] = q[i]; |
| 138 | } |
| 139 | for (j = 0; j < n; ++j) // 0..n-1 |
| 140 | { |
| 141 | Matrix::Row UDj(UD,j); |
| 142 | Matrix::const_Row Gj(G,j); |
| 143 | for (i = 0; i < Nq; ++i) // 0..Nq-1 |
| 144 | UDj[i+n] = Gj[i]; |
| 145 | } |
| 146 | |
| 147 | // U=Fx*U and diagonals retrieved |
| 148 | for (j = n-1; j > 0; --j) // n-1..1 |
| 149 | { |
| 150 | // Prepare d(0)..d(j) as temporary |
| 151 | for (i = 0; i <= j; ++i) // 0..j |
| 152 | d[i] = Float(UD(i,j)); // ISSUE mixed type proxy assignment |
| 153 | |
| 154 | // Lower triangle of UD is implicitly empty |
| 155 | for (i = 0; i < n; ++i) // 0..n-1 |
| 156 | { |
| 157 | Matrix::Row UDi(UD,i); |
| 158 | Matrix::const_Row Fxi(Fx,i); |
| 159 | UDi[j] = Fxi[j]; |
| 160 | for (k = 0; k < j; ++k) // 0..j-1 |
| 161 | UDi[j] += Fxi[k] * d[k]; |
| 162 | } |
| 163 | } |
| 164 | d[0] = Float(UD(0,0)); // ISSUE mixed type proxy assignment |
| 165 | |
| 166 | // Complete U = Fx*U |
| 167 | for (j = 0; j < n; ++j) // 0..n-1 |
nothing calls this directly
no test coverage detected