OMCompiler/Compiler/NSimCode/NSimJacobian.mo
| Line | Branch | Exec | Source |
|---|---|---|---|
| 1 | /* | ||
| 2 | * This file is part of OpenModelica. | ||
| 3 | * | ||
| 4 | * Copyright (c) 1998-2026, Open Source Modelica Consortium (OSMC), | ||
| 5 | * c/o Linköpings universitet, Department of Computer and Information Science, | ||
| 6 | * SE-58183 Linköping, Sweden. | ||
| 7 | * | ||
| 8 | * All rights reserved. | ||
| 9 | * | ||
| 10 | * THIS PROGRAM IS PROVIDED UNDER THE TERMS OF AGPL VERSION 3 LICENSE OR | ||
| 11 | * THIS OSMC PUBLIC LICENSE (OSMC-PL) VERSION 1.8. | ||
| 12 | * ANY USE, REPRODUCTION OR DISTRIBUTION OF THIS PROGRAM CONSTITUTES | ||
| 13 | * RECIPIENT'S ACCEPTANCE OF THE OSMC PUBLIC LICENSE OR THE GNU AGPL | ||
| 14 | * VERSION 3, ACCORDING TO RECIPIENTS CHOICE. | ||
| 15 | * | ||
| 16 | * The OpenModelica software and the OSMC (Open Source Modelica Consortium) | ||
| 17 | * Public License (OSMC-PL) are obtained from OSMC, either from the above | ||
| 18 | * address, from the URLs: | ||
| 19 | * http://www.openmodelica.org or | ||
| 20 | * https://github.com/OpenModelica/ or | ||
| 21 | * http://www.ida.liu.se/projects/OpenModelica, | ||
| 22 | * and in the OpenModelica distribution. | ||
| 23 | * | ||
| 24 | * GNU AGPL version 3 is obtained from: | ||
| 25 | * https://www.gnu.org/licenses/licenses.html#GPL | ||
| 26 | * | ||
| 27 | * This program is distributed WITHOUT ANY WARRANTY; without | ||
| 28 | * even the implied warranty of MERCHANTABILITY or FITNESS | ||
| 29 | * FOR A PARTICULAR PURPOSE, EXCEPT AS EXPRESSLY SET FORTH | ||
| 30 | * IN THE BY RECIPIENT SELECTED SUBSIDIARY LICENSE CONDITIONS OF OSMC-PL. | ||
| 31 | * | ||
| 32 | * See the full OSMC Public License conditions for more details. | ||
| 33 | * | ||
| 34 | */ | ||
| 35 | |||
| 36 | encapsulated package NSimJacobian | ||
| 37 | "file: NSimJacobian.mo | ||
| 38 | package: NSimJacobian | ||
| 39 | description: This file contains the functions for creating simcode jaobians and sparsity patterns. | ||
| 40 | " | ||
| 41 | |||
| 42 | public | ||
| 43 | // NF imports | ||
| 44 | import ComponentRef = NFComponentRef; | ||
| 45 | import Dimension = NFDimension; | ||
| 46 | import Expression = NFExpression; | ||
| 47 | import NFInstNode.InstNode; | ||
| 48 | import Subscript = NFSubscript; | ||
| 49 | import Type = NFType; | ||
| 50 | |||
| 51 | // Backend imports | ||
| 52 | import Adjacency = NBAdjacency; | ||
| 53 | import NBAdjacency.Dependency; | ||
| 54 | import BackendDAE = NBackendDAE; | ||
| 55 | import NBEquation.{Equation, Iterator, EquationPointer, EquationPointers, EqData}; | ||
| 56 | import BEquation = NBEquation; | ||
| 57 | import NBVariable.{VariablePointers, VarData}; | ||
| 58 | import BVariable = NBVariable; | ||
| 59 | import Variable = NFVariable; | ||
| 60 | import Jacobian = NBJacobian; | ||
| 61 | import Partition = NBPartition; | ||
| 62 | |||
| 63 | // SimCode imports | ||
| 64 | import SimCodeUtil = NSimCodeUtil; | ||
| 65 | import SimCode = NSimCode; | ||
| 66 | import SimGenericCall = NSimGenericCall; | ||
| 67 | import NSimCode.Identifier; | ||
| 68 | import SimStrongComponent = NSimStrongComponent; | ||
| 69 | import NSimVar.{SimVar, SimVars, VarType, ConvertMemo}; | ||
| 70 | |||
| 71 | // Old SimCode imports | ||
| 72 | import OldSimCode = SimCode; | ||
| 73 | import OldBackendDAE = BackendDAE; | ||
| 74 | |||
| 75 | // Util imports | ||
| 76 | import StringUtil; | ||
| 77 | |||
| 78 | uniontype SparsityRow | ||
| 79 | record SPARSITY_ROW | ||
| 80 | ComponentRef equation_name "only for debugging"; | ||
| 81 | list<SimGenericCall.SimIterator> equation_iterators; | ||
| 82 | list<tuple<ComponentRef, Dependency, Boolean /*true=repeated*/>> dependencies; | ||
| 83 | list<ComponentRef> solved_crefs; | ||
| 84 | end SPARSITY_ROW; | ||
| 85 | |||
| 86 | function create | ||
| 87 | input ComponentRef equation_name; | ||
| 88 | input Iterator equation_iterator; | ||
| 89 | input UnorderedMap<ComponentRef, Dependency> dependencies; | ||
| 90 | input UnorderedSet<ComponentRef> repetitions; | ||
| 91 | input list<ComponentRef> solved_crefs; | ||
| 92 | output SparsityRow row; | ||
| 93 | protected | ||
| 94 | list<ComponentRef> crefs; | ||
| 95 | list<Dependency> deps; | ||
| 96 | list<Boolean> reps; | ||
| 97 | algorithm | ||
| 98 | 871 | crefs := UnorderedMap.keyList(dependencies); | |
| 99 | 871 | deps := UnorderedMap.valueList(dependencies); | |
| 100 |
6/6✓ Branch 0 taken 2115 times.
✓ Branch 1 taken 871 times.
✓ Branch 2 taken 2115 times.
✓ Branch 3 taken 871 times.
✓ Branch 5 taken 2053 times.
✓ Branch 6 taken 62 times.
|
2986 | reps := list(UnorderedSet.contains(cref, repetitions) for cref in crefs); |
| 101 | |||
| 102 | // add whole subscripts for code gen purposes | ||
| 103 |
4/4✓ Branch 0 taken 2115 times.
✓ Branch 1 taken 871 times.
✓ Branch 2 taken 2115 times.
✓ Branch 3 taken 871 times.
|
2986 | crefs := list(ComponentRef.fillSubscripts(cref) for cref in crefs); |
| 104 | |||
| 105 |
13/14✓ Branch 0 taken 2115 times.
✓ Branch 1 taken 871 times.
✓ Branch 2 taken 2115 times.
✓ Branch 3 taken 871 times.
✓ Branch 4 taken 2115 times.
✓ Branch 5 taken 871 times.
✓ Branch 6 taken 2115 times.
✓ Branch 7 taken 871 times.
✗ Branch 9 not taken.
✓ Branch 10 taken 871 times.
✓ Branch 12 taken 873 times.
✓ Branch 13 taken 871 times.
✓ Branch 14 taken 873 times.
✓ Branch 15 taken 871 times.
|
3859 | row := SPARSITY_ROW( |
| 106 | equation_name = equation_name, | ||
| 107 | equation_iterators = SimGenericCall.SimIterator.fromIterator(equation_iterator), | ||
| 108 | dependencies = list((cref, dep, rep) threaded for cref in crefs, dep in deps, rep in reps), | ||
| 109 | solved_crefs = list(ComponentRef.fillSubscripts(cref) for cref in solved_crefs)); | ||
| 110 | end create; | ||
| 111 | |||
| 112 | function convert | ||
| 113 | input SparsityRow row; | ||
| 114 | output OldSimCode.SparsityRow oldrow; | ||
| 115 | algorithm | ||
| 116 |
12/12✓ Branch 0 taken 101 times.
✓ Branch 1 taken 2416 times.
✓ Branch 2 taken 101 times.
✓ Branch 3 taken 2416 times.
✓ Branch 5 taken 6121 times.
✓ Branch 6 taken 2416 times.
✓ Branch 7 taken 6121 times.
✓ Branch 8 taken 2416 times.
✓ Branch 15 taken 2417 times.
✓ Branch 16 taken 2416 times.
✓ Branch 17 taken 2417 times.
✓ Branch 18 taken 2416 times.
|
11055 | oldrow := OldSimCode.SPARSITY_ROW( |
| 117 | equation_name = ComponentRef.toDAE(row.equation_name), | ||
| 118 | equation_iterators = list(SimGenericCall.SimIterator.convert(iter) for iter in row.equation_iterators), | ||
| 119 | dependencies = list((ComponentRef.toDAE(Util.tuple31(tpl)), Dependency.convert(Util.tuple32(tpl)), Util.tuple33(tpl)) for tpl in row.dependencies), | ||
| 120 | solved_crefs = list(ComponentRef.toDAE(cref) for cref in row.solved_crefs) | ||
| 121 | ); | ||
| 122 | end convert; | ||
| 123 | |||
| 124 | function toString | ||
| 125 | input SparsityRow row; | ||
| 126 | output String str; | ||
| 127 | function dependencyString | ||
| 128 | input tuple<ComponentRef, Dependency, Boolean> tpl; | ||
| 129 | output String str = "(" + ComponentRef.toString(Util.tuple31(tpl)) + ", " + Dependency.toString(Util.tuple32(tpl)) + ", " + boolString(Util.tuple33(tpl)) + ")"; | ||
| 130 | end dependencyString; | ||
| 131 | algorithm | ||
| 132 | 76 | str := ComponentRef.toString(row.equation_name) + " ... " + List.toString(row.solved_crefs, ComponentRef.toString) + " ... " + List.toString(row.dependencies, dependencyString); | |
| 133 | end toString; | ||
| 134 | |||
| 135 | function mergeDuplicateRows | ||
| 136 | "Index reduction can occasionally produce more residual equations than | ||
| 137 | there are physical Jacobian rows (numberOfResultVars, the runtime's | ||
| 138 | fixed row capacity) -- e.g. static_IR.mos, where an over-determined | ||
| 139 | algebraic subsystem yields several equivalent constraints that all end | ||
| 140 | up tagged with the same solved_crefs. In that situation resizable | ||
| 141 | Jacobian codegen (CodegenC.tpl) would still assign each extra equation | ||
| 142 | its own incrementing row index, writing past the runtime's | ||
| 143 | numberOfResultVars-sized buffers. | ||
| 144 | This is NOT the same as a legitimately coupled multi-row system (e.g. | ||
| 145 | algebraicLoop.mos), where several genuinely different equations validly | ||
| 146 | share the same solved_crefs set (one row each, all coupled to the same | ||
| 147 | variables) with row count already matching numberOfResultVars -- such | ||
| 148 | cases must be left untouched, hence the row-count guard below." | ||
| 149 | input list<SparsityRow> rows_in; | ||
| 150 | input Integer numberOfResultVars; | ||
| 151 | output list<SparsityRow> rows_out; | ||
| 152 | protected | ||
| 153 | UnorderedMap<String, SparsityRow> row_map; | ||
| 154 | list<String> order; | ||
| 155 | String key; | ||
| 156 | SparsityRow existing, merged; | ||
| 157 | algorithm | ||
| 158 |
1/2✗ Branch 1 not taken.
✓ Branch 2 taken 167 times.
|
167 | if listLength(rows_in) <= numberOfResultVars then |
| 159 | rows_out := rows_in; | ||
| 160 | else | ||
| 161 | ✗ | row_map := UnorderedMap.new<SparsityRow>(stringHashDjb2, stringEq); | |
| 162 | order := {}; | ||
| 163 | ✗ | for row in rows_in loop | |
| 164 | ✗ | key := stringDelimitList(list(ComponentRef.toString(c) for c in row.solved_crefs), ","); | |
| 165 | ✗ | if UnorderedMap.contains(key, row_map) then | |
| 166 | ✗ | existing := UnorderedMap.getSafe(key, row_map, sourceInfo()); | |
| 167 | merged := existing; | ||
| 168 | ✗ | merged.dependencies := List.unionOnTrue(existing.dependencies, row.dependencies, dependencyCrefEqual); | |
| 169 | ✗ | UnorderedMap.add(key, merged, row_map); | |
| 170 | else | ||
| 171 | ✗ | UnorderedMap.add(key, row, row_map); | |
| 172 | order := key :: order; | ||
| 173 | end if; | ||
| 174 | end for; | ||
| 175 | ✗ | order := listReverse(order); | |
| 176 | ✗ | rows_out := list(UnorderedMap.getSafe(k, row_map, sourceInfo()) for k in order); | |
| 177 | end if; | ||
| 178 | end mergeDuplicateRows; | ||
| 179 | |||
| 180 | function mergeScalarRows | ||
| 181 | "Merges rows that solve the same scalar variables. The rows of these variables | ||
| 182 | are their absolute positions, so every equation of an algebraic loop, which | ||
| 183 | solves all its unknowns, would add the same rows again with its dependencies. | ||
| 184 | This multiplied the non-zeros (with duplicates) and the generated code by the | ||
| 185 | number of loop equations. The merged row has the union of the dependencies. | ||
| 186 | Rows with equation iterators or array solved variables keep their own rows." | ||
| 187 | input list<SparsityRow> rows_in; | ||
| 188 | output list<SparsityRow> rows_out = {}; | ||
| 189 | protected | ||
| 190 | UnorderedMap<String, SparsityRow> row_map = UnorderedMap.new<SparsityRow>(stringHashDjb2, stringEq); | ||
| 191 | type DepSet = UnorderedSet<ComponentRef>; | ||
| 192 | UnorderedMap<String, DepSet> dep_map = UnorderedMap.new<DepSet>(stringHashDjb2, stringEq); | ||
| 193 | UnorderedSet<ComponentRef> deps; | ||
| 194 | String key; | ||
| 195 | SparsityRow merged; | ||
| 196 | algorithm | ||
| 197 |
2/2✓ Branch 0 taken 871 times.
✓ Branch 1 taken 167 times.
|
1038 | for row in rows_in loop |
| 198 |
2/2✓ Branch 1 taken 796 times.
✓ Branch 2 taken 75 times.
|
871 | if isScalarRow(row) then |
| 199 |
4/4✓ Branch 0 taken 798 times.
✓ Branch 1 taken 796 times.
✓ Branch 2 taken 798 times.
✓ Branch 3 taken 796 times.
|
1594 | key := stringDelimitList(list(ComponentRef.toString(c) for c in row.solved_crefs), ","); |
| 200 |
2/2✓ Branch 1 taken 1 time.
✓ Branch 2 taken 795 times.
|
796 | if UnorderedMap.contains(key, row_map) then |
| 201 | 1 | merged := UnorderedMap.getSafe(key, row_map, sourceInfo()); | |
| 202 | 1 | deps := UnorderedMap.getSafe(key, dep_map, sourceInfo()); | |
| 203 |
2/2✓ Branch 0 taken 2 times.
✓ Branch 1 taken 1 time.
|
3 | for dep in row.dependencies loop |
| 204 |
1/2✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
|
2 | if not UnorderedSet.contains(Util.tuple31(dep), deps) then |
| 205 | ✗ | UnorderedSet.add(Util.tuple31(dep), deps); | |
| 206 | ✗ | merged.dependencies := dep :: merged.dependencies; | |
| 207 | end if; | ||
| 208 | end for; | ||
| 209 | 1 | UnorderedMap.add(key, merged, row_map); | |
| 210 | else | ||
| 211 | 795 | UnorderedMap.add(key, row, row_map); | |
| 212 |
4/4✓ Branch 0 taken 1815 times.
✓ Branch 1 taken 795 times.
✓ Branch 2 taken 1815 times.
✓ Branch 3 taken 795 times.
|
2610 | UnorderedMap.add(key, UnorderedSet.fromList(list(Util.tuple31(dep) for dep in row.dependencies), ComponentRef.hash, ComponentRef.isEqual), dep_map); |
| 213 | rows_out := row :: rows_out; | ||
| 214 | end if; | ||
| 215 | else | ||
| 216 | rows_out := row :: rows_out; | ||
| 217 | end if; | ||
| 218 | end for; | ||
| 219 | // replace the first occurrences by the merged rows | ||
| 220 |
10/10✓ Branch 0 taken 870 times.
✓ Branch 1 taken 167 times.
✓ Branch 2 taken 870 times.
✓ Branch 3 taken 167 times.
✓ Branch 5 taken 795 times.
✓ Branch 6 taken 75 times.
✓ Branch 7 taken 796 times.
✓ Branch 8 taken 795 times.
✓ Branch 9 taken 796 times.
✓ Branch 10 taken 795 times.
|
1833 | rows_out := listReverse(list(if isScalarRow(row) then UnorderedMap.getSafe(stringDelimitList(list(ComponentRef.toString(c) for c in row.solved_crefs), ","), row_map, sourceInfo()) else row for row in rows_out)); |
| 221 | end mergeScalarRows; | ||
| 222 | |||
| 223 | function isScalarRow | ||
| 224 | "true if the rows of the solved variables are absolute positions: no equation | ||
| 225 | iterators and only index subscripts" | ||
| 226 | input SparsityRow row; | ||
| 227 | output Boolean b; | ||
| 228 | algorithm | ||
| 229 |
5/6✓ Branch 0 taken 1645 times.
✓ Branch 1 taken 96 times.
✓ Branch 2 taken 1645 times.
✗ Branch 3 not taken.
✓ Branch 5 taken 54 times.
✓ Branch 6 taken 1591 times.
|
1741 | b := listEmpty(row.equation_iterators) and not listEmpty(row.solved_crefs) |
| 230 | and List.all(row.solved_crefs, crefHasOnlyIndexSubscripts); | ||
| 231 | end isScalarRow; | ||
| 232 | |||
| 233 | function crefHasOnlyIndexSubscripts | ||
| 234 | input ComponentRef cref; | ||
| 235 | output Boolean b = List.all(ComponentRef.subscriptsAllFlat(cref), Subscript.isIndex); | ||
| 236 | end crefHasOnlyIndexSubscripts; | ||
| 237 | |||
| 238 | function sortByResultVars | ||
| 239 | "Orders the rows like the result variables. The runtime reads row i of the | ||
| 240 | pattern as result variable i, but the rows come in equation order. | ||
| 241 | Leaves the rows unchanged if a solved variable is not a result variable." | ||
| 242 | input list<SparsityRow> rows_in; | ||
| 243 | input list<SimVar> resVars; | ||
| 244 | output list<SparsityRow> rows_out = rows_in; | ||
| 245 | protected | ||
| 246 | UnorderedMap<ComponentRef, Integer> first_index = UnorderedMap.new<Integer>(ComponentRef.hash, ComponentRef.isEqual); | ||
| 247 | ComponentRef name; | ||
| 248 | list<tuple<Integer, Integer, SparsityRow>> keyed = {}; | ||
| 249 | Integer key, pos = 0; | ||
| 250 | algorithm | ||
| 251 |
2/2✓ Branch 0 taken 1441 times.
✓ Branch 1 taken 164 times.
|
1605 | for sv in resVars loop |
| 252 | 1441 | name := ComponentRef.stripSubscriptsAll(sv.name); | |
| 253 | 1441 | UnorderedMap.add(name, intMin(sv.index, UnorderedMap.getOrDefault(name, first_index, sv.index)), first_index); | |
| 254 | end for; | ||
| 255 | |||
| 256 |
2/2✓ Branch 0 taken 859 times.
✓ Branch 1 taken 164 times.
|
1023 | for row in rows_in loop |
| 257 | key := -1; | ||
| 258 |
2/2✓ Branch 1 taken 860 times.
✓ Branch 2 taken 859 times.
|
1719 | for cref in row.solved_crefs loop |
| 259 | () := match UnorderedMap.get(ComponentRef.stripSubscriptsAll(cref), first_index) | ||
| 260 | local | ||
| 261 | Integer idx; | ||
| 262 | case SOME(idx) algorithm | ||
| 263 |
2/2✓ Branch 0 taken 1 time.
✓ Branch 1 taken 859 times.
|
860 | key := if key < 0 then idx else intMin(key, idx); |
| 264 | then (); | ||
| 265 | else (); | ||
| 266 | end match; | ||
| 267 | end for; | ||
| 268 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 859 times.
|
859 | if key < 0 then |
| 269 | ✗ | return; | |
| 270 | end if; | ||
| 271 | 859 | keyed := (key, pos, row) :: keyed; | |
| 272 | 859 | pos := pos + 1; | |
| 273 | end for; | ||
| 274 | |||
| 275 | 164 | keyed := List.sort(keyed, rowKeyGreater); | |
| 276 |
4/4✓ Branch 0 taken 859 times.
✓ Branch 1 taken 164 times.
✓ Branch 2 taken 859 times.
✓ Branch 3 taken 164 times.
|
1023 | rows_out := list(Util.tuple33(t) for t in keyed); |
| 277 | end sortByResultVars; | ||
| 278 | |||
| 279 | function rowKeyGreater | ||
| 280 | input tuple<Integer, Integer, SparsityRow> row1; | ||
| 281 | input tuple<Integer, Integer, SparsityRow> row2; | ||
| 282 | output Boolean b; | ||
| 283 | protected | ||
| 284 | Integer key1, key2, pos1, pos2; | ||
| 285 | algorithm | ||
| 286 | 2714 | (key1, pos1, _) := row1; | |
| 287 | 2714 | (key2, pos2, _) := row2; | |
| 288 |
3/4✓ Branch 0 taken 2695 times.
✓ Branch 1 taken 19 times.
✓ Branch 2 taken 2695 times.
✗ Branch 3 not taken.
|
2714 | b := key1 > key2 or (key1 == key2 and pos1 > pos2); |
| 289 | end rowKeyGreater; | ||
| 290 | |||
| 291 | function dependencyCrefEqual | ||
| 292 | input tuple<ComponentRef, Dependency, Boolean> dep1; | ||
| 293 | input tuple<ComponentRef, Dependency, Boolean> dep2; | ||
| 294 | output Boolean b; | ||
| 295 | algorithm | ||
| 296 | ✗ | b := ComponentRef.isEqual(Util.tuple31(dep1), Util.tuple31(dep2)); | |
| 297 | end dependencyCrefEqual; | ||
| 298 | end SparsityRow; | ||
| 299 | |||
| 300 | uniontype Sparsity | ||
| 301 | record SPARSITY | ||
| 302 | list<SparsityRow> rows; | ||
| 303 | end SPARSITY; | ||
| 304 | |||
| 305 | record EMPTY | ||
| 306 | end EMPTY; | ||
| 307 | |||
| 308 | function create | ||
| 309 | input Adjacency.Matrix mat; | ||
| 310 | input list<SimVar> resVars; | ||
| 311 | input Boolean isAdjoint; | ||
| 312 | output Sparsity sparsity; | ||
| 313 | protected | ||
| 314 | list<SparsityRow> rows; | ||
| 315 | algorithm | ||
| 316 | sparsity := match mat | ||
| 317 | case Adjacency.SPARSITY() algorithm | ||
| 318 |
13/14✓ Branch 0 taken 871 times.
✓ Branch 1 taken 167 times.
✓ Branch 3 taken 871 times.
✓ Branch 4 taken 167 times.
✓ Branch 6 taken 871 times.
✓ Branch 7 taken 167 times.
✓ Branch 9 taken 871 times.
✓ Branch 10 taken 167 times.
✓ Branch 12 taken 871 times.
✓ Branch 13 taken 167 times.
✓ Branch 15 taken 871 times.
✓ Branch 16 taken 167 times.
✗ Branch 18 not taken.
✓ Branch 19 taken 167 times.
|
6228 | rows := SparsityRow.mergeDuplicateRows( |
| 319 | list(SparsityRow.create(e, i, d, r, s) threaded for e in mat.equation_names, i in mat.equation_iterators, d in mat.dependencies, r in mat.repetitions, s in mat.solved_crefs), | ||
| 320 | SimVars.numScalarElems(resVars)); | ||
| 321 | 167 | rows := SparsityRow.mergeScalarRows(rows); | |
| 322 |
2/2✓ Branch 0 taken 164 times.
✓ Branch 1 taken 3 times.
|
167 | if not isAdjoint then |
| 323 | 164 | rows := SparsityRow.sortByResultVars(rows, resVars); | |
| 324 | end if; | ||
| 325 | 167 | then SPARSITY(rows); | |
| 326 | case Adjacency.EMPTY() then EMPTY(); | ||
| 327 | |||
| 328 | else algorithm | ||
| 329 | ✗ | Error.addMessage(Error.INTERNAL_ERROR,{getInstanceName() + " can only handle sparsity or empty matrices but got:\n" + Adjacency.Matrix.toString(mat)}); | |
| 330 | ✗ | then fail(); | |
| 331 | end match; | ||
| 332 | end create; | ||
| 333 | |||
| 334 | function convert | ||
| 335 | input Sparsity sparsity; | ||
| 336 | output OldSimCode.Sparsity oldsparsity; | ||
| 337 | algorithm | ||
| 338 | oldsparsity := match sparsity | ||
| 339 |
4/4✓ Branch 0 taken 2416 times.
✓ Branch 1 taken 369 times.
✓ Branch 2 taken 2416 times.
✓ Branch 3 taken 369 times.
|
2785 | case SPARSITY() then OldSimCode.SPARSITY(list(SparsityRow.convert(row) for row in sparsity.rows)); |
| 340 | case EMPTY() then OldSimCode.EMPTY(); | ||
| 341 | end match; | ||
| 342 | end convert; | ||
| 343 | |||
| 344 | function toString | ||
| 345 | input Sparsity sparsity; | ||
| 346 | output String str = StringUtil.headline_3("Resizable Sparsity Pattern"); | ||
| 347 | algorithm | ||
| 348 | str := match sparsity | ||
| 349 | 14 | case SPARSITY() then str + List.toString(sparsity.rows, SparsityRow.toString, List.Style.NEWLINE) + "\n"; | |
| 350 | ✗ | else str + " -- EMPTY -- \n"; | |
| 351 | end match; | ||
| 352 | end toString; | ||
| 353 | end Sparsity; | ||
| 354 | |||
| 355 | uniontype SimJacobian | ||
| 356 | record SIM_JAC | ||
| 357 | String name "unique matrix name"; | ||
| 358 | Integer jacobianIndex "unique jacobian index"; | ||
| 359 | Integer partitionIndex "index of partition it belongs to"; | ||
| 360 | Integer numberOfResultVars "corresponds to the number of rows"; | ||
| 361 | list<SimStrongComponent.Block> columnEqns "column equations equals in size to column vars"; | ||
| 362 | list<SimStrongComponent.Block> constantEqns "List of constant equations independent of seed variables"; | ||
| 363 | list<SimVar> columnVars "all column vars, none results vars index -1, the other corresponding to rows index"; | ||
| 364 | list<SimVar> seedVars "corresponds to the number of columns"; | ||
| 365 | Sparsity sparsityMatrix "new sparsity pattern"; | ||
| 366 | list<SimGenericCall> generic_loop_calls "Generic for-loop and array calls"; | ||
| 367 | Option<UnorderedMap<ComponentRef, SimVar>> jac_map "hash table for cref -> simVar"; | ||
| 368 | Boolean isAdjoint "indicates if this is an adjoint jacobian"; | ||
| 369 | Boolean isBidirectional "indicates if this jacobian is part of a bidirectional pair"; | ||
| 370 | Integer adjointJacobianIndex "index of the adjoint jacobian for bidirectional (-1 if not bidirectional)"; | ||
| 371 | String adjointMatrixName "matrix name of the adjoint jacobian for bidirectional"; | ||
| 372 | end SIM_JAC; | ||
| 373 | |||
| 374 | function toString | ||
| 375 | input SimJacobian simJac; | ||
| 376 | output String str = ""; | ||
| 377 | algorithm | ||
| 378 | str := match simJac | ||
| 379 | case SIM_JAC() algorithm | ||
| 380 |
2/2✓ Branch 1 taken 58 times.
✓ Branch 2 taken 14 times.
|
72 | if isEmpty(simJac) then |
| 381 | 58 | str := StringUtil.headline_2("[EMPTY] SimCode Jacobian " + simJac.name + "(idx = " + intString(simJac.jacobianIndex) + ", partition = " + intString(simJac.partitionIndex) + ")") + "\n"; | |
| 382 | else | ||
| 383 | 14 | str := StringUtil.headline_2("SimCode Jacobian " + simJac.name + "(idx = " + intString(simJac.jacobianIndex) + ", partition = " + intString(simJac.jacobianIndex) + ")") + "\n"; | |
| 384 | 14 | str := str + StringUtil.headline_4("SeedVars (size = " + intString(listLength(simJac.seedVars)) + ")"); | |
| 385 |
2/2✓ Branch 0 taken 79 times.
✓ Branch 1 taken 14 times.
|
93 | for var in simJac.seedVars loop |
| 386 | 79 | str := str + SimVar.toString(var, " ") + "\n"; | |
| 387 | end for; | ||
| 388 | 14 | str := str + "\n" + StringUtil.headline_4("TmpVars (size = " + intString(listLength(simJac.columnVars)) + ")"); | |
| 389 |
2/2✓ Branch 0 taken 26 times.
✓ Branch 1 taken 14 times.
|
40 | for var in simJac.columnVars loop |
| 390 | 26 | str := str + SimVar.toString(var, " ") + "\n"; | |
| 391 | end for; | ||
| 392 | // TODO: print list of ResultVars | ||
| 393 | 14 | str := str + "\n" + StringUtil.headline_4("ResultVars (size = " + intString(simJac.numberOfResultVars) + ")"); | |
| 394 | |||
| 395 | // TODO: count equations properly, e.g. linear systems are falsely counted as a single equation | ||
| 396 | 14 | str := str + "\n" + StringUtil.headline_3("Column Equations (size = " + intString(listLength(simJac.columnEqns)) + ")"); | |
| 397 |
2/2✓ Branch 0 taken 89 times.
✓ Branch 1 taken 14 times.
|
103 | for eq in simJac.columnEqns loop |
| 398 | 89 | str := str + SimStrongComponent.Block.toString(eq, " "); | |
| 399 | end for; | ||
| 400 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 14 times.
|
14 | if not listEmpty(simJac.constantEqns) then |
| 401 | ✗ | str := str + StringUtil.headline_3("Constant Equations"); | |
| 402 | ✗ | for eq in simJac.constantEqns loop | |
| 403 | ✗ | str := str + SimStrongComponent.Block.toString(eq, " "); | |
| 404 | end for; | ||
| 405 | end if; | ||
| 406 | 14 | str := str + "\n" + Sparsity.toString(simJac.sparsityMatrix); | |
| 407 | |||
| 408 |
2/2✓ Branch 0 taken 1 time.
✓ Branch 1 taken 13 times.
|
14 | if not listEmpty(simJac.generic_loop_calls) then |
| 409 | 1 | str := str + StringUtil.headline_3("Generic Calls"); | |
| 410 | 1 | str := str + List.toString(simJac.generic_loop_calls, SimGenericCall.toString, List.Style.NEWLINE_INDENT); | |
| 411 | end if; | ||
| 412 | 14 | str := str + "\n"; | |
| 413 | end if; | ||
| 414 | then str; | ||
| 415 | else algorithm | ||
| 416 | ✗ | Error.addMessage(Error.INTERNAL_ERROR, {getInstanceName() + " failed."}); | |
| 417 | ✗ | then fail(); | |
| 418 | end match; | ||
| 419 | end toString; | ||
| 420 | |||
| 421 | function isEmpty | ||
| 422 | input SimJacobian simJac; | ||
| 423 | output Boolean b; | ||
| 424 | algorithm | ||
| 425 | b := match simJac | ||
| 426 | 72 | case SIM_JAC() then simJac.numberOfResultVars == 0; | |
| 427 | else false; | ||
| 428 | end match; | ||
| 429 | end isEmpty; | ||
| 430 | |||
| 431 | /* | ||
| 432 | function fromSystems | ||
| 433 | input list<System.System> partitions; | ||
| 434 | output Option<SimJacobian> simJacobian; | ||
| 435 | input output SimCode.SimCodeIndices indices; | ||
| 436 | protected | ||
| 437 | list<BackendDAE> jacobians = {}; | ||
| 438 | algorithm | ||
| 439 | for partition in partitions loop | ||
| 440 | if isSome(partition.jacobian) then | ||
| 441 | jacobians := Util.getOption(partition.jacobian) :: jacobians; | ||
| 442 | end if; | ||
| 443 | end for; | ||
| 444 | |||
| 445 | if listEmpty(jacobians) then | ||
| 446 | simJacobian := NONE(); | ||
| 447 | else | ||
| 448 | (simJacobian, indices) := create(Jacobian.combine(jacobians, "A"), indices); | ||
| 449 | end if; | ||
| 450 | end fromSystems; | ||
| 451 | |||
| 452 | function fromSystemsSparsity | ||
| 453 | input list<System.System> partitions; | ||
| 454 | input output Option<SimJacobian> simJacobian; | ||
| 455 | input UnorderedMap<ComponentRef, SimVar> sim_map; | ||
| 456 | input output SimCode.SimCodeIndices indices; | ||
| 457 | algorithm | ||
| 458 | (simJacobian, indices) := match (partitions, simJacobian) | ||
| 459 | local | ||
| 460 | BackendDAE jacobian; | ||
| 461 | |||
| 462 | case (_, NONE()) then (NONE(), indices); | ||
| 463 | case ({System.SYSTEM(jacobian = NONE())}, _) then (NONE(), indices); | ||
| 464 | case ({System.SYSTEM(jacobian = SOME(jacobian))}, _) then createSparsity(jacobian, Util.getOption(simJacobian), sim_map, indices); | ||
| 465 | else algorithm | ||
| 466 | Error.addMessage(Error.INTERNAL_ERROR,{getInstanceName() + " failed! Partitioned partitions are not yet supported by this function."}); | ||
| 467 | then fail(); | ||
| 468 | |||
| 469 | end match; | ||
| 470 | end fromSystemsSparsity; | ||
| 471 | */ | ||
| 472 | function create | ||
| 473 | input BackendDAE jacobian; | ||
| 474 | output Option<SimJacobian> simJacobian; | ||
| 475 | input output SimCode.SimCodeIndices indices; | ||
| 476 | input UnorderedMap<ComponentRef, SimVar> simcode_map; | ||
| 477 | algorithm | ||
| 478 | simJacobian := match jacobian | ||
| 479 | local | ||
| 480 | // dummy map for strong component creation (no alias possible here) | ||
| 481 | UnorderedMap<ComponentRef, SimVar> dummy_sim_map = UnorderedMap.new<SimVar>(ComponentRef.hash, ComponentRef.isEqual); | ||
| 482 | UnorderedMap<ComponentRef, SimStrongComponent.Block> dummy_eqn_map = UnorderedMap.new<SimStrongComponent.Block>(ComponentRef.hash, ComponentRef.isEqual); | ||
| 483 | SimStrongComponent.Block columnEqn; | ||
| 484 | list<SimStrongComponent.Block> columnEqns = {}; | ||
| 485 | VarData varData; | ||
| 486 | list<Pointer<Variable>> seed_lst, res_lst, tmp_lst; | ||
| 487 | list<SimVar> seedVars, resVars, tmpVars; | ||
| 488 | UnorderedMap<ComponentRef, SimVar> jac_map; | ||
| 489 | SimJacobian jac; | ||
| 490 | UnorderedMap<Identifier, Integer> sim_map; | ||
| 491 | list<SimGenericCall> generic_loop_calls; | ||
| 492 | UnorderedMap<ComponentRef, Integer> min_sub_map; | ||
| 493 | UnorderedMap<ComponentRef, SimVar> min_sv_map; | ||
| 494 | |||
| 495 | case BackendDAE.JACOBIAN(varData = varData as BVariable.VAR_DATA_JAC()) algorithm | ||
| 496 | // temporarily save the generic call map from simcode to recover it afterwards | ||
| 497 | // we use a local map to have seperated generic call lists for each jacobian | ||
| 498 | 167 | sim_map := indices.generic_call_map; | |
| 499 | 167 | indices.generic_call_map := UnorderedMap.new<Integer>(Identifier.hash, Identifier.isEqual); | |
| 500 |
3/4✗ Branch 0 not taken.
✓ Branch 1 taken 167 times.
✓ Branch 2 taken 88 times.
✓ Branch 3 taken 79 times.
|
1482 | for i in arrayLength(jacobian.comps):-1:1 loop |
| 501 | 1148 | (columnEqn, indices, _) := SimStrongComponent.Block.fromStrongComponent(jacobian.comps[i], indices, NBPartition.Kind.JAC, dummy_sim_map, dummy_eqn_map); | |
| 502 | columnEqns := columnEqn :: columnEqns; | ||
| 503 | end for; | ||
| 504 | |||
| 505 | // extract generic loop calls and put the old generic call map back | ||
| 506 |
4/4✓ Branch 1 taken 62 times.
✓ Branch 2 taken 167 times.
✓ Branch 3 taken 62 times.
✓ Branch 4 taken 167 times.
|
229 | generic_loop_calls := list(SimGenericCall.fromIdentifier(tpl) for tpl in UnorderedMap.toList(indices.generic_call_map)); |
| 507 | 167 | indices.generic_call_map := sim_map; | |
| 508 | |||
| 509 |
2/2✓ Branch 1 taken 126 times.
✓ Branch 2 taken 41 times.
|
167 | if Flags.getConfigBool(Flags.SIM_CODE_SCALARIZE) then |
| 510 | 126 | seed_lst := VariablePointers.toList(VariablePointers.scalarize(varData.seedVars)); | |
| 511 | 126 | res_lst := VariablePointers.toList(VariablePointers.scalarize(varData.resultVars)); | |
| 512 | 126 | tmp_lst := VariablePointers.toList(VariablePointers.scalarize(varData.tmpVars)); | |
| 513 | else | ||
| 514 | 41 | seed_lst := VariablePointers.toList(varData.seedVars); | |
| 515 | 41 | res_lst := VariablePointers.toList(varData.resultVars); | |
| 516 | 41 | tmp_lst := VariablePointers.toList(varData.tmpVars); | |
| 517 | end if; | ||
| 518 | |||
| 519 | // the runtime reads column and row i of the ODE Jacobian as state i | ||
| 520 |
4/4✓ Branch 0 taken 81 times.
✓ Branch 1 taken 86 times.
✓ Branch 2 taken 78 times.
✓ Branch 3 taken 3 times.
|
167 | if jacobian.jacType == NBJacobian.JacobianType.ODE and not jacobian.isAdjoint then |
| 521 | 78 | seed_lst := sortByStateIndex(seed_lst, simcode_map); | |
| 522 | 78 | res_lst := sortByStateIndex(res_lst, simcode_map); | |
| 523 | end if; | ||
| 524 | |||
| 525 | // column and seed var indices always start at 0. Without scalarization the | ||
| 526 | // Jacobian has no index map like the simulation variables, so the index of an | ||
| 527 | // array variable is the position of its first element in seedVars, resultVars | ||
| 528 | // and tmpVars, and in the columns and rows of the sparsity pattern | ||
| 529 | 167 | seedVars := SimVar.createList(seed_lst, VarType.SIMULATION, NSimCode.EMPTY_SIM_CODE_INDICES(), not Flags.getConfigBool(Flags.SIM_CODE_SCALARIZE)); | |
| 530 | 167 | resVars := SimVar.createList(res_lst, VarType.SIMULATION, NSimCode.EMPTY_SIM_CODE_INDICES(), not Flags.getConfigBool(Flags.SIM_CODE_SCALARIZE)); | |
| 531 | 167 | tmpVars := SimVar.createList(tmp_lst, VarType.SIMULATION, NSimCode.EMPTY_SIM_CODE_INDICES(), not Flags.getConfigBool(Flags.SIM_CODE_SCALARIZE)); | |
| 532 | |||
| 533 | 167 | jac_map := UnorderedMap.new<SimVar>(ComponentRef.hash, ComponentRef.isEqual, listLength(seedVars) + listLength(resVars) + listLength(tmpVars)); | |
| 534 | 167 | SimCodeUtil.addListSimCodeMap(seedVars, jac_map); | |
| 535 | 167 | SimCodeUtil.addListSimCodeMap(resVars, jac_map); | |
| 536 | 167 | SimCodeUtil.addListSimCodeMap(tmpVars, jac_map); | |
| 537 | |||
| 538 | // For per-element seed groups that start above index [1] (partial-slice | ||
| 539 | // NLS iter vars), add a virtual SimVar keyed on the base cref (no | ||
| 540 | // subscripts). Templates look up $SEED.x when they see $SEED.x[$i1] | ||
| 541 | // with an iterator subscript; the virtual SimVar's index encodes the | ||
| 542 | // offset so (&seedVars[index])[$i1-1] reaches the correct element. | ||
| 543 | // | ||
| 544 | // Two-pass: first find the MINIMUM outer subscript per base-cref group | ||
| 545 | // (the hash-ordered seedVars traversal does not guarantee ascending | ||
| 546 | // module order), then create virtual SimVars using that minimum seed so | ||
| 547 | // the base index is always non-negative. | ||
| 548 | 167 | min_sub_map := UnorderedMap.new<Integer>(ComponentRef.hash, ComponentRef.isEqual); | |
| 549 | 167 | min_sv_map := UnorderedMap.new<SimVar>(ComponentRef.hash, ComponentRef.isEqual); | |
| 550 | |||
| 551 |
2/2✓ Branch 0 taken 1451 times.
✓ Branch 1 taken 167 times.
|
1618 | for sv in seedVars loop |
| 552 |
2/2✓ Branch 1 taken 775 times.
✓ Branch 2 taken 676 times.
|
1451 | if ComponentRef.hasSubscripts(sv.name) then |
| 553 | () := match ComponentRef.outermostIntegerSubscript(sv.name) | ||
| 554 | local | ||
| 555 | ComponentRef base_cref; | ||
| 556 | Integer sub_val, cur_min; | ||
| 557 | case sub_val guard sub_val > 1 | ||
| 558 | algorithm | ||
| 559 | 636 | base_cref := ComponentRef.stripSubscriptsAll(sv.name); | |
| 560 |
2/2✓ Branch 1 taken 563 times.
✓ Branch 2 taken 73 times.
|
636 | if UnorderedMap.contains(base_cref, min_sub_map) then |
| 561 | 563 | cur_min := UnorderedMap.getSafe(base_cref, min_sub_map, sourceInfo()); | |
| 562 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 563 times.
|
563 | if sub_val < cur_min then |
| 563 | ✗ | UnorderedMap.add(base_cref, sub_val, min_sub_map); | |
| 564 | ✗ | UnorderedMap.add(base_cref, sv, min_sv_map); | |
| 565 | end if; | ||
| 566 | else | ||
| 567 | 73 | UnorderedMap.add(base_cref, sub_val, min_sub_map); | |
| 568 | 73 | UnorderedMap.add(base_cref, sv, min_sv_map); | |
| 569 | end if; | ||
| 570 | then (); | ||
| 571 | else (); | ||
| 572 | end match; | ||
| 573 | end if; | ||
| 574 | end for; | ||
| 575 | |||
| 576 |
2/2✓ Branch 1 taken 73 times.
✓ Branch 2 taken 167 times.
|
240 | for tpl in UnorderedMap.toList(min_sv_map) loop |
| 577 | () := match tpl | ||
| 578 | local | ||
| 579 | ComponentRef base_cref; | ||
| 580 | SimVar min_sv, virtual_sv; | ||
| 581 | Integer cur_min; | ||
| 582 | case (base_cref, min_sv) | ||
| 583 | algorithm | ||
| 584 | 73 | cur_min := UnorderedMap.getSafe(base_cref, min_sub_map, sourceInfo()); | |
| 585 |
1/2✓ Branch 1 taken 73 times.
✗ Branch 2 not taken.
|
73 | if not UnorderedMap.contains(base_cref, jac_map) then |
| 586 | virtual_sv := min_sv; | ||
| 587 | 73 | virtual_sv.name := base_cref; | |
| 588 | 73 | virtual_sv.index := min_sv.index - crefFlatOffset(min_sv.name); | |
| 589 | 73 | virtual_sv.arrayCref := NONE(); | |
| 590 | 73 | UnorderedMap.add(base_cref, virtual_sv, jac_map); | |
| 591 | end if; | ||
| 592 | then (); | ||
| 593 | end match; | ||
| 594 | end for; | ||
| 595 | |||
| 596 | 167 | jac := SIM_JAC( | |
| 597 | name = jacobian.name, | ||
| 598 | jacobianIndex = indices.jacobianIndex, | ||
| 599 | partitionIndex = 0, | ||
| 600 | numberOfResultVars = SimVars.numScalarElems(resVars), | ||
| 601 | columnEqns = columnEqns, | ||
| 602 | constantEqns = {}, | ||
| 603 | columnVars = tmpVars, | ||
| 604 | seedVars = seedVars, | ||
| 605 | sparsityMatrix = Sparsity.create(jacobian.sparsity, resVars, jacobian.isAdjoint), | ||
| 606 | generic_loop_calls = generic_loop_calls, | ||
| 607 | jac_map = SOME(jac_map), | ||
| 608 | isAdjoint = jacobian.isAdjoint, | ||
| 609 | isBidirectional = false, | ||
| 610 | adjointJacobianIndex = -1, | ||
| 611 | adjointMatrixName = "" | ||
| 612 | ); | ||
| 613 | |||
| 614 | 167 | indices.jacobianIndex := indices.jacobianIndex + 1; | |
| 615 | simJacobian := SOME(jac); | ||
| 616 | then simJacobian; | ||
| 617 | |||
| 618 | else algorithm | ||
| 619 | ✗ | Error.addMessage(Error.INTERNAL_ERROR, {getInstanceName() + " failed."}); | |
| 620 | ✗ | then fail(); | |
| 621 | end match; | ||
| 622 | end create; | ||
| 623 | |||
| 624 | function createSimulationJacobian | ||
| 625 | input list<Partition.Partition> partitions; | ||
| 626 | output SimJacobian simJac; | ||
| 627 | output SimJacobian simJacAdjoint; | ||
| 628 | input output SimCode.SimCodeIndices simCodeIndices; | ||
| 629 | input UnorderedMap<ComponentRef, SimVar> simcode_map; | ||
| 630 | protected | ||
| 631 | list<BackendDAE> jacobians = {}, jacobiansAdjoint = {}; | ||
| 632 | BackendDAE simJacobian, simJacobianAdjoint; | ||
| 633 | Option<SimJacobian> simJac_opt, simJacAdj_opt; | ||
| 634 | Option<BackendDAE> jacobian, jacobianAdjoint; | ||
| 635 | algorithm | ||
| 636 |
2/2✓ Branch 0 taken 83 times.
✓ Branch 1 taken 188 times.
|
271 | for partition in partitions loop |
| 637 | // save jacobian if existent | ||
| 638 | 83 | jacobian := Partition.Partition.getJacobian(partition); | |
| 639 |
3/4✗ Branch 0 not taken.
✓ Branch 1 taken 83 times.
✓ Branch 2 taken 80 times.
✓ Branch 3 taken 3 times.
|
83 | if isSome(jacobian) then |
| 640 | 80 | jacobians := Util.getOption(jacobian) :: jacobians; | |
| 641 | end if; | ||
| 642 | 83 | jacobianAdjoint := Partition.Partition.getJacobianAdjoint(partition); | |
| 643 |
3/4✗ Branch 0 not taken.
✓ Branch 1 taken 83 times.
✓ Branch 2 taken 4 times.
✓ Branch 3 taken 79 times.
|
83 | if isSome(jacobianAdjoint) then |
| 644 | 4 | jacobiansAdjoint := Util.getOption(jacobianAdjoint) :: jacobiansAdjoint; | |
| 645 | end if; | ||
| 646 | end for; | ||
| 647 | |||
| 648 | // create empty jacobian as fallback | ||
| 649 | // ToDo: handle simCodeIndices correctly here | ||
| 650 |
2/2✓ Branch 0 taken 109 times.
✓ Branch 1 taken 79 times.
|
188 | if listEmpty(jacobians) then |
| 651 | 109 | (simJac, simCodeIndices) := SimJacobian.empty("A", simCodeIndices); | |
| 652 | else | ||
| 653 | 79 | simJacobian := Jacobian.combine(jacobians, "A"); | |
| 654 | 79 | (simJac_opt, simCodeIndices) := SimJacobian.create(simJacobian, simCodeIndices, simcode_map); | |
| 655 |
2/4✗ Branch 0 not taken.
✓ Branch 1 taken 79 times.
✓ Branch 2 taken 79 times.
✗ Branch 3 not taken.
|
79 | if isSome(simJac_opt) then |
| 656 | 79 | simJac := Util.getOption(simJac_opt); | |
| 657 | else | ||
| 658 | ✗ | (simJac, simCodeIndices) := SimJacobian.empty("A", simCodeIndices); | |
| 659 | end if; | ||
| 660 | end if; | ||
| 661 | |||
| 662 | // create empty adjoint jacobian as fallback | ||
| 663 |
2/2✓ Branch 0 taken 185 times.
✓ Branch 1 taken 3 times.
|
188 | if listEmpty(jacobiansAdjoint) then |
| 664 | 185 | (simJacAdjoint, simCodeIndices) := SimJacobian.empty("ADJ", simCodeIndices); | |
| 665 | else | ||
| 666 | 3 | simJacobianAdjoint := Jacobian.combine(jacobiansAdjoint, "ADJ"); | |
| 667 | 3 | (simJacAdj_opt, simCodeIndices) := SimJacobian.create(simJacobianAdjoint, simCodeIndices, simcode_map); | |
| 668 |
2/4✗ Branch 0 not taken.
✓ Branch 1 taken 3 times.
✓ Branch 2 taken 3 times.
✗ Branch 3 not taken.
|
3 | if isSome(simJacAdj_opt) then |
| 669 | 3 | simJacAdjoint := Util.getOption(simJacAdj_opt); | |
| 670 | else | ||
| 671 | ✗ | (simJacAdjoint, simCodeIndices) := SimJacobian.empty("ADJ", simCodeIndices); | |
| 672 | end if; | ||
| 673 | end if; | ||
| 674 | |||
| 675 | // Link forward and adjoint for bidirectional mode | ||
| 676 |
3/4✓ Branch 1 taken 1 time.
✓ Branch 2 taken 187 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 1 time.
|
188 | if Flags.getConfigString(Flags.GENERATE_DYNAMIC_JACOBIAN) == "bidirectional" then |
| 677 | simJac := match simJac | ||
| 678 | case SIM_JAC() algorithm | ||
| 679 | 1 | simJac.isBidirectional := true; | |
| 680 | 1 | simJac.adjointJacobianIndex := match simJacAdjoint case SIM_JAC() then simJacAdjoint.jacobianIndex; else -1; end match; | |
| 681 | 1 | simJac.adjointMatrixName := match simJacAdjoint case SIM_JAC() then simJacAdjoint.name; else ""; end match; | |
| 682 | then simJac; | ||
| 683 | else simJac; | ||
| 684 | end match; | ||
| 685 | end if; | ||
| 686 | end createSimulationJacobian; | ||
| 687 | |||
| 688 | function createOptimizationJacobian | ||
| 689 | input list<Partition.Partition> partitions; | ||
| 690 | output SimJacobian simJacLfg; | ||
| 691 | output SimJacobian simJacMrf; | ||
| 692 | output SimJacobian simJacR0; | ||
| 693 | input output SimCode.SimCodeIndices simCodeIndices; | ||
| 694 | input UnorderedMap<ComponentRef, SimVar> simcode_map; | ||
| 695 | protected | ||
| 696 | list<BackendDAE> jacobiansLfg = {}, jacobiansMrf = {}, jacobiansR0 = {}; | ||
| 697 | BackendDAE simJacobianLfg, simJacobianMrf, simJacobianR0; | ||
| 698 | Option<SimJacobian> simJacLfg_opt, simJacMrf_opt, simJacR0_opt; | ||
| 699 | Option<BackendDAE> jacobianLfg, jacobianMrf, jacobianR0; | ||
| 700 | algorithm | ||
| 701 |
2/2✓ Branch 0 taken 83 times.
✓ Branch 1 taken 188 times.
|
271 | for partition in partitions loop |
| 702 | // collect | ||
| 703 | 83 | jacobianLfg := Partition.Partition.getJacobianLfg(partition); | |
| 704 |
2/4✗ Branch 0 not taken.
✓ Branch 1 taken 83 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 83 times.
|
83 | if isSome(jacobianLfg) then |
| 705 | ✗ | jacobiansLfg := Util.getOption(jacobianLfg) :: jacobiansLfg; | |
| 706 | end if; | ||
| 707 | 83 | jacobianMrf := Partition.Partition.getJacobianMrf(partition); | |
| 708 |
2/4✗ Branch 0 not taken.
✓ Branch 1 taken 83 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 83 times.
|
83 | if isSome(jacobianMrf) then |
| 709 | ✗ | jacobiansMrf := Util.getOption(jacobianMrf) :: jacobiansMrf; | |
| 710 | end if; | ||
| 711 | 83 | jacobianR0 := Partition.Partition.getJacobianR0(partition); | |
| 712 |
2/4✗ Branch 0 not taken.
✓ Branch 1 taken 83 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 83 times.
|
83 | if isSome(jacobianR0) then |
| 713 | ✗ | jacobiansR0 := Util.getOption(jacobianR0) :: jacobiansR0; | |
| 714 | end if; | ||
| 715 | end for; | ||
| 716 | |||
| 717 | // create empty Lfg jacobian as fallback | ||
| 718 |
1/2✓ Branch 0 taken 188 times.
✗ Branch 1 not taken.
|
188 | if listEmpty(jacobiansLfg) then |
| 719 | 188 | (simJacLfg, simCodeIndices) := SimJacobian.empty("OPT_LFG", simCodeIndices); | |
| 720 | else | ||
| 721 | ✗ | simJacobianLfg := Jacobian.combine(jacobiansLfg, "OPT_LFG"); | |
| 722 | ✗ | (simJacLfg_opt, simCodeIndices) := SimJacobian.create(simJacobianLfg, simCodeIndices, simcode_map); | |
| 723 | ✗ | if isSome(simJacLfg_opt) then | |
| 724 | ✗ | simJacLfg := Util.getOption(simJacLfg_opt); | |
| 725 | else | ||
| 726 | ✗ | (simJacLfg, simCodeIndices) := SimJacobian.empty("OPT_LFG", simCodeIndices); | |
| 727 | end if; | ||
| 728 | end if; | ||
| 729 | |||
| 730 | // create empty Mrf jacobian as fallback | ||
| 731 |
1/2✓ Branch 0 taken 188 times.
✗ Branch 1 not taken.
|
188 | if listEmpty(jacobiansMrf) then |
| 732 | 188 | (simJacMrf, simCodeIndices) := SimJacobian.empty("OPT_MRF", simCodeIndices); | |
| 733 | else | ||
| 734 | ✗ | simJacobianMrf := Jacobian.combine(jacobiansMrf, "OPT_MRF"); | |
| 735 | ✗ | (simJacMrf_opt, simCodeIndices) := SimJacobian.create(simJacobianMrf, simCodeIndices, simcode_map); | |
| 736 | ✗ | if isSome(simJacMrf_opt) then | |
| 737 | ✗ | simJacMrf := Util.getOption(simJacMrf_opt); | |
| 738 | else | ||
| 739 | ✗ | (simJacMrf, simCodeIndices) := SimJacobian.empty("OPT_MRF", simCodeIndices); | |
| 740 | end if; | ||
| 741 | end if; | ||
| 742 | |||
| 743 | // create empty R0 jacobian as fallback | ||
| 744 |
1/2✓ Branch 0 taken 188 times.
✗ Branch 1 not taken.
|
188 | if listEmpty(jacobiansR0) then |
| 745 | 188 | (simJacR0, simCodeIndices) := SimJacobian.empty("OPT_R0", simCodeIndices); | |
| 746 | else | ||
| 747 | ✗ | simJacobianR0 := Jacobian.combine(jacobiansR0, "OPT_R0"); | |
| 748 | ✗ | (simJacR0_opt, simCodeIndices) := SimJacobian.create(simJacobianR0, simCodeIndices, simcode_map); | |
| 749 | ✗ | if isSome(simJacR0_opt) then | |
| 750 | ✗ | simJacR0 := Util.getOption(simJacR0_opt); | |
| 751 | else | ||
| 752 | ✗ | (simJacR0, simCodeIndices) := SimJacobian.empty("OPT_R0", simCodeIndices); | |
| 753 | end if; | ||
| 754 | end if; | ||
| 755 | end createOptimizationJacobian; | ||
| 756 | |||
| 757 | function empty | ||
| 758 | input String name = ""; | ||
| 759 | output SimJacobian emptyJac = EMPTY_SIM_JAC; | ||
| 760 | input output SimCode.SimCodeIndices indices; | ||
| 761 | algorithm | ||
| 762 | emptyJac := match emptyJac | ||
| 763 | case SIM_JAC() algorithm | ||
| 764 | 1798 | emptyJac.name := name; | |
| 765 | 1798 | emptyJac.jacobianIndex := indices.jacobianIndex; | |
| 766 |
1/2✓ Branch 0 taken 1798 times.
✗ Branch 1 not taken.
|
1798 | indices.jacobianIndex := indices.jacobianIndex + 1; |
| 767 | then emptyJac; | ||
| 768 | else algorithm | ||
| 769 | ✗ | Error.addMessage(Error.INTERNAL_ERROR, {getInstanceName() + " failed."}); | |
| 770 | ✗ | then fail(); | |
| 771 | end match; | ||
| 772 | end empty; | ||
| 773 | |||
| 774 | function getJacobianBlocks | ||
| 775 | input SimJacobian jacobian; | ||
| 776 | output list<SimStrongComponent.Block> blcks; | ||
| 777 | algorithm | ||
| 778 | blcks := match jacobian | ||
| 779 | 1880 | case SIM_JAC() then listAppend(jacobian.constantEqns, jacobian.columnEqns); | |
| 780 | else algorithm | ||
| 781 | ✗ | Error.addMessage(Error.INTERNAL_ERROR, {getInstanceName() + " failed."}); | |
| 782 | ✗ | then fail(); | |
| 783 | end match; | ||
| 784 | end getJacobianBlocks; | ||
| 785 | |||
| 786 | function getJacobiansBlocks | ||
| 787 | input list<SimJacobian> jacobians; | ||
| 788 | output list<SimStrongComponent.Block> blcks = {}; | ||
| 789 | algorithm | ||
| 790 |
2/2✓ Branch 0 taken 1880 times.
✓ Branch 1 taken 188 times.
|
2068 | for jacobian in jacobians loop |
| 791 | 1880 | blcks := listAppend(getJacobianBlocks(jacobian), blcks); | |
| 792 | end for; | ||
| 793 | end getJacobiansBlocks; | ||
| 794 | |||
| 795 | function getJacobianHT | ||
| 796 | input SimJacobian jacobian; | ||
| 797 | output Option<UnorderedMap<ComponentRef, SimVar>> jac_map; | ||
| 798 | algorithm | ||
| 799 | jac_map := match jacobian | ||
| 800 | ✗ | case SIM_JAC() then jacobian.jac_map; | |
| 801 | else algorithm | ||
| 802 | ✗ | Error.addMessage(Error.INTERNAL_ERROR, {getInstanceName() + " failed."}); | |
| 803 | ✗ | then fail(); | |
| 804 | end match; | ||
| 805 | end getJacobianHT; | ||
| 806 | |||
| 807 | function convert | ||
| 808 | input SimJacobian simJac; | ||
| 809 | output OldSimCode.JacobianMatrix oldJac; | ||
| 810 | protected | ||
| 811 | OldSimCode.JacobianColumn oldJacCol; | ||
| 812 | algorithm | ||
| 813 | oldJac := match simJac | ||
| 814 | local | ||
| 815 | ConvertMemo memo; | ||
| 816 | case SIM_JAC() algorithm | ||
| 817 | 2167 | memo := SimVar.newConvertMemo(listLength(simJac.seedVars) + listLength(simJac.columnVars)); | |
| 818 |
6/8✓ Branch 0 taken 3707 times.
✓ Branch 1 taken 2167 times.
✓ Branch 2 taken 3707 times.
✓ Branch 3 taken 2167 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2167 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 2167 times.
|
5874 | oldJacCol := OldSimCode.JAC_COLUMN( |
| 819 | columnEqns = list(SimStrongComponent.Block.convert(blck) for blck in simJac.columnEqns), | ||
| 820 | columnVars = SimVar.convertListMemo(simJac.columnVars, memo), | ||
| 821 | numberOfResultVars = simJac.numberOfResultVars, | ||
| 822 | constantEqns = list(SimStrongComponent.Block.convert(blck) for blck in simJac.constantEqns) | ||
| 823 | ); | ||
| 824 | |||
| 825 |
4/4✓ Branch 0 taken 212 times.
✓ Branch 1 taken 2167 times.
✓ Branch 2 taken 212 times.
✓ Branch 3 taken 2167 times.
|
2379 | oldJac := OldSimCode.JAC_MATRIX( |
| 826 | columns = {oldJacCol}, | ||
| 827 | seedVars = SimVar.convertListMemo(simJac.seedVars, memo), | ||
| 828 | matrixName = simJac.name, | ||
| 829 | sparsityMatrix = Sparsity.convert(simJac.sparsityMatrix), | ||
| 830 | sparsity = {}, | ||
| 831 | sparsityT = {}, | ||
| 832 | nonlinear = {}, | ||
| 833 | nonlinearT = {}, | ||
| 834 | coloredCols = {}, | ||
| 835 | coloredRows = {}, | ||
| 836 | maxColorCols = 0, | ||
| 837 | jacobianIndex = simJac.jacobianIndex, | ||
| 838 | partitionIndex = simJac.partitionIndex, | ||
| 839 | generic_loop_calls = list(SimGenericCall.convert(gc) for gc in simJac.generic_loop_calls), | ||
| 840 | crefsHT = Util.applyOption(simJac.jac_map, function SimCodeUtil.convertSimCodeMap(memo = memo)), | ||
| 841 | isAdjoint = simJac.isAdjoint, | ||
| 842 | isBidirectional = simJac.isBidirectional, | ||
| 843 | adjointJacobianIndex = simJac.adjointJacobianIndex, | ||
| 844 | adjointMatrixName = simJac.adjointMatrixName | ||
| 845 | ); | ||
| 846 | then oldJac; | ||
| 847 | |||
| 848 | else algorithm | ||
| 849 | ✗ | Error.addMessage(Error.INTERNAL_ERROR, {getInstanceName() + " failed."}); | |
| 850 | ✗ | then fail(); | |
| 851 | end match; | ||
| 852 | end convert; | ||
| 853 | end SimJacobian; | ||
| 854 | |||
| 855 | constant SimJacobian EMPTY_SIM_JAC = SIM_JAC("", 0, 0, 0, {}, {}, {}, {}, | ||
| 856 | Sparsity.EMPTY(), {}, NONE(), false, false, -1, ""); | ||
| 857 | |||
| 858 | protected function collectNodeSubDimPairsOuterFirst | ||
| 859 | "Collect (sub_val, dim_size) pairs for a single CREF node's subscripts and | ||
| 860 | dimension list, in outer→inner (natural subscript) order. Only INTEGER | ||
| 861 | subscripts are included; non-integer subscripts consume their dimension slot | ||
| 862 | without emitting a pair. Always has an else case." | ||
| 863 | input list<Subscript> subs; | ||
| 864 | input list<Dimension> dims; | ||
| 865 | output list<tuple<Integer, Integer>> pairs; | ||
| 866 | algorithm | ||
| 867 | pairs := match (subs, dims) | ||
| 868 | local | ||
| 869 | Subscript s; | ||
| 870 | list<Subscript> rest_subs; | ||
| 871 | Dimension d; | ||
| 872 | list<Dimension> rest_dims; | ||
| 873 | Integer v, v3; | ||
| 874 | list<tuple<Integer, Integer>> rest_pairs; | ||
| 875 | case ({}, _) then {}; | ||
| 876 | 20 | case ({Subscript.INDEX(index = Expression.INTEGER(v3))}, {}) then {(v3, 1)}; | |
| 877 | // InstNode.getType(node) can come back without array dims for a seed's own | ||
| 878 | // node (e.g. an INIT-partition per-element seed cref, whose leaf node's | ||
| 879 | // declared type resolves to a bare scalar even though it indexes an array | ||
| 880 | // -- see the ODE partition's equivalent cref, whose node keeps its proper | ||
| 881 | // Real[n] type, for the same conceptual variable). When this is the node's | ||
| 882 | // ONLY subscript, falling all the way through to `{}` (as the general | ||
| 883 | // multi-subscript case below still does) would silently drop this | ||
| 884 | // dimension's contribution to the flat offset entirely -- exactly the | ||
| 885 | // "sizeCols/size mismatch" class of bug this whole file exists to avoid. | ||
| 886 | // A dim_size of 1 is always safe here: this pair is the innermost/only | ||
| 887 | // one from this node, so nothing multiplies by it going outward, and | ||
| 888 | // (v3-1)*1 is exactly the correct offset contribution. | ||
| 889 | case (_, {}) then {}; | ||
| 890 | case (s :: rest_subs, d :: rest_dims) | ||
| 891 | algorithm | ||
| 892 | 66 | rest_pairs := collectNodeSubDimPairsOuterFirst(rest_subs, rest_dims); | |
| 893 | then | ||
| 894 | match s | ||
| 895 | local Integer v2; | ||
| 896 | 66 | case Subscript.INDEX(index = Expression.INTEGER(v2)) then (v2, Dimension.size(d)) :: rest_pairs; | |
| 897 | else rest_pairs; | ||
| 898 | end match; | ||
| 899 | else {}; | ||
| 900 | end match; | ||
| 901 | end collectNodeSubDimPairsOuterFirst; | ||
| 902 | |||
| 903 | protected function crefSubDimPairsLeafToRoot | ||
| 904 | "Collect (integer_subscript_value, dimension_size) pairs from the innermost | ||
| 905 | (leaf) CREF node to the outermost (root), in inner-first order globally. | ||
| 906 | Within each node, multiple subscripts are handled in inner→outer order | ||
| 907 | (reversed from the outer→inner subscript list so the accumulation works). | ||
| 908 | Uses InstNode.getType(node) (the pre-subscript type) rather than cref.ty | ||
| 909 | (the post-subscript element type) so that record-field subscripts like | ||
| 910 | module[2] (where cref.ty = Module, not Module[10]) recover the full array | ||
| 911 | dimension needed to compute correct flat offsets." | ||
| 912 | input ComponentRef cref; | ||
| 913 | output list<tuple<Integer, Integer>> pairs; | ||
| 914 | algorithm | ||
| 915 | pairs := match cref | ||
| 916 | local | ||
| 917 | list<Subscript> subs; | ||
| 918 | Type node_ty; | ||
| 919 | ComponentRef rest; | ||
| 920 | list<tuple<Integer, Integer>> rest_pairs, node_pairs; | ||
| 921 | list<Dimension> dims; | ||
| 922 | Integer own; | ||
| 923 | case ComponentRef.CREF(subscripts = subs, restCref = rest) | ||
| 924 | algorithm | ||
| 925 | 204 | rest_pairs := crefSubDimPairsLeafToRoot(rest); | |
| 926 | // Use the node's own declared type (before applying these subscripts) | ||
| 927 | // so record-valued fields also expose their array dimensions. | ||
| 928 | 204 | node_ty := InstNode.getType(ComponentRef.node(cref)); | |
| 929 | // the type of a record field can have the dimensions of its parents lifted in front | ||
| 930 | 204 | dims := Type.arrayDims(node_ty); | |
| 931 | 204 | own := listLength(dims) - listLength(ComponentRef.subscriptsAllFlat(rest)); | |
| 932 |
4/4✓ Branch 1 taken 181 times.
✓ Branch 2 taken 23 times.
✓ Branch 4 taken 14 times.
✓ Branch 5 taken 167 times.
|
204 | if own >= listLength(subs) and own < listLength(dims) then |
| 933 | 14 | dims := List.lastN(dims, own); | |
| 934 | end if; | ||
| 935 | // Collect this node's pairs outer-first, then reverse to get inner-first. | ||
| 936 | // Append rest_pairs (which are from the outer/restCref direction) after. | ||
| 937 | 204 | node_pairs := listReverse(collectNodeSubDimPairsOuterFirst(subs, dims)); | |
| 938 | 204 | then | |
| 939 | listAppend(node_pairs, rest_pairs); | ||
| 940 | else {}; | ||
| 941 | end match; | ||
| 942 | end crefSubDimPairsLeafToRoot; | ||
| 943 | |||
| 944 | protected function sortByStateIndex | ||
| 945 | "Orders seed or result variables like the simulation variables they belong to. | ||
| 946 | Unchanged if one of them is not a simulation variable." | ||
| 947 | input list<Pointer<Variable>> vars; | ||
| 948 | input UnorderedMap<ComponentRef, SimVar> simcode_map; | ||
| 949 | output list<Pointer<Variable>> sorted = vars; | ||
| 950 | protected | ||
| 951 | list<tuple<Integer, Pointer<Variable>>> keyed = {}; | ||
| 952 | Option<SimVar> osv; | ||
| 953 | SimVar sv; | ||
| 954 | algorithm | ||
| 955 |
2/2✓ Branch 0 taken 1180 times.
✓ Branch 1 taken 156 times.
|
1336 | for v in vars loop |
| 956 | 1180 | osv := UnorderedMap.get(stripRootNode(BVariable.getVarName(v)), simcode_map); | |
| 957 |
2/4✗ Branch 0 not taken.
✓ Branch 1 taken 1180 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 1180 times.
|
1180 | if isNone(osv) then |
| 958 | ✗ | return; | |
| 959 | end if; | ||
| 960 | 1180 | SOME(sv) := osv; | |
| 961 | 1180 | keyed := (sv.index, v) :: keyed; | |
| 962 | end for; | ||
| 963 | 156 | keyed := List.sort(keyed, indexGreater); | |
| 964 |
4/4✓ Branch 0 taken 1180 times.
✓ Branch 1 taken 156 times.
✓ Branch 2 taken 1180 times.
✓ Branch 3 taken 156 times.
|
1336 | sorted := list(Util.tuple22(t) for t in keyed); |
| 965 | end sortByStateIndex; | ||
| 966 | |||
| 967 | protected function indexGreater | ||
| 968 | input tuple<Integer, Pointer<Variable>> t1; | ||
| 969 | input tuple<Integer, Pointer<Variable>> t2; | ||
| 970 | output Boolean b = Util.tuple21(t1) > Util.tuple21(t2); | ||
| 971 | end indexGreater; | ||
| 972 | |||
| 973 | protected function stripRootNode | ||
| 974 | "$SEED_ODE_JAC.x -> x" | ||
| 975 | input ComponentRef cref; | ||
| 976 | output ComponentRef stripped; | ||
| 977 | algorithm | ||
| 978 | stripped := match cref | ||
| 979 | case ComponentRef.CREF(restCref = ComponentRef.EMPTY()) then ComponentRef.EMPTY(); | ||
| 980 | case ComponentRef.CREF() algorithm | ||
| 981 | 2084 | cref.restCref := stripRootNode(cref.restCref); | |
| 982 | then cref; | ||
| 983 | else cref; | ||
| 984 | end match; | ||
| 985 | end stripRootNode; | ||
| 986 | |||
| 987 | protected function crefFlatOffset | ||
| 988 | "Compute the 0-based flat array index encoded by all integer subscripts in | ||
| 989 | the cref chain, using row-major (C-style) layout. For example, x[2][1] in | ||
| 990 | a 10×10 array returns (2-1)*10 + (1-1) = 10. Used to compute the virtual | ||
| 991 | SimVar index so that (&seedVars[virtual.idx])[(i1-1)*d2+...+(iN-1)] maps | ||
| 992 | correctly to the actual seed index for every valid subscript combination." | ||
| 993 | input ComponentRef cref; | ||
| 994 | output Integer offset = 0; | ||
| 995 | protected | ||
| 996 | list<tuple<Integer, Integer>> pairs; | ||
| 997 | Integer inner_prod = 1, sub_val, dim_sz; | ||
| 998 | algorithm | ||
| 999 | // pairs is in inner-first order: process each pair with accumulated inner_prod. | ||
| 1000 | // flat = sum_i (sub[i]-1) * product_of_dims_more_inner_than_i | ||
| 1001 | 73 | pairs := crefSubDimPairsLeafToRoot(cref); | |
| 1002 |
2/2✓ Branch 0 taken 86 times.
✓ Branch 1 taken 73 times.
|
159 | for pair in pairs loop |
| 1003 | 86 | (sub_val, dim_sz) := pair; | |
| 1004 | 86 | offset := offset + (sub_val - 1) * inner_prod; | |
| 1005 | 86 | inner_prod := inner_prod * dim_sz; | |
| 1006 | end for; | ||
| 1007 | end crefFlatOffset; | ||
| 1008 | |||
| 1009 | annotation(__OpenModelica_Interface="nbackend"); | ||
| 1010 | end NSimJacobian; | ||
| 1011 |