OMCompiler/Compiler/NBackEnd/Modules/3_Post/NBJacobian.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 NBJacobian | ||
| 37 | "file: NBJacobian.mo | ||
| 38 | package: NBJacobian | ||
| 39 | description: This file contains the functions to create and manipulate jacobians. | ||
| 40 | The main type is inherited from NBackendDAE.mo | ||
| 41 | NOTE: There is no real jacobian type, it is a BackendDAE. | ||
| 42 | " | ||
| 43 | |||
| 44 | public | ||
| 45 | import BackendDAE = NBackendDAE; | ||
| 46 | import Module = NBModule; | ||
| 47 | import NBEquation; | ||
| 48 | import NBVariable; | ||
| 49 | |||
| 50 | protected | ||
| 51 | // OF imports | ||
| 52 | import Absyn.Path; | ||
| 53 | import DAE; | ||
| 54 | |||
| 55 | // NF imports | ||
| 56 | import ComponentRef = NFComponentRef; | ||
| 57 | import Algorithm = NFAlgorithm; | ||
| 58 | import Expression = NFExpression; | ||
| 59 | import NFFunction.Function; | ||
| 60 | import Statement = NFStatement; | ||
| 61 | import Subscript = NFSubscript; | ||
| 62 | import List; | ||
| 63 | import Operator = NFOperator; | ||
| 64 | import SimplifyExp = NFSimplifyExp; | ||
| 65 | import Type = NFType; | ||
| 66 | import Variable = NFVariable; | ||
| 67 | import Scalarize = NFScalarize; | ||
| 68 | import NFInstNode.InstNode; | ||
| 69 | |||
| 70 | // Backend imports | ||
| 71 | import NFBackendExtension.BackendInfo; | ||
| 72 | import Adjacency = NBAdjacency; | ||
| 73 | import NBAdjacency.Mapping; | ||
| 74 | import BEquation = NBEquation; | ||
| 75 | import BVariable = NBVariable; | ||
| 76 | import Differentiate = NBDifferentiate; | ||
| 77 | import NBDifferentiate.{DifferentiationArguments, DifferentiationType}; | ||
| 78 | import NBEquation.{Equation, Iterator, EquationPointers, EqData}; | ||
| 79 | import Jacobian = NBackendDAE.BackendDAE; | ||
| 80 | import Matching = NBMatching; | ||
| 81 | import Partition = NBPartition; | ||
| 82 | import Replacements = NBReplacements; | ||
| 83 | import Slice = NBSlice; | ||
| 84 | import Sorting = NBSorting; | ||
| 85 | import StrongComponent = NBStrongComponent; | ||
| 86 | import Tearing = NBTearing; | ||
| 87 | import NFOperator.{MathClassification, SizeClassification}; | ||
| 88 | import NBVariable.{VariablePointer, VariablePointers, VarData}; | ||
| 89 | |||
| 90 | // Sparsity-pattern graph coloring, shared with the old backend. | ||
| 91 | import Coloring; | ||
| 92 | |||
| 93 | // Util imports | ||
| 94 | import StringUtil; | ||
| 95 | import PointerWeak; | ||
| 96 | import UnorderedMap; | ||
| 97 | import UnorderedSet; | ||
| 98 | import Util; | ||
| 99 | |||
| 100 | public | ||
| 101 | type JacobianType = enumeration(ODE, DAE, LS, NLS, OPT_LFG, OPT_MRF, OPT_R0); | ||
| 102 | |||
| 103 | function isDynamic | ||
| 104 | "is the jacobian used for integration (-> true) | ||
| 105 | or solving algebraic systems (-> false)?" | ||
| 106 | input JacobianType jacType; | ||
| 107 | output Boolean b; | ||
| 108 | algorithm | ||
| 109 | b := match jacType | ||
| 110 | case JacobianType.ODE then true; | ||
| 111 | case JacobianType.DAE then true; | ||
| 112 | case JacobianType.OPT_LFG then true; | ||
| 113 | case JacobianType.OPT_MRF then true; | ||
| 114 | case JacobianType.OPT_R0 then true; | ||
| 115 | else false; | ||
| 116 | end match; | ||
| 117 | end isDynamic; | ||
| 118 | |||
| 119 | function main | ||
| 120 | "Wrapper function for any jacobian function. This will be called during | ||
| 121 | simulation and gets the corresponding subfunction from Config." | ||
| 122 | extends Module.wrapper; | ||
| 123 | input Partition.Kind kind; | ||
| 124 | protected | ||
| 125 | constant Module.jacobianInterface func = getModule(); | ||
| 126 | algorithm | ||
| 127 | bdae := match bdae | ||
| 128 | local | ||
| 129 | String name "Context name for jacobian"; | ||
| 130 | VariablePointers knowns "Variable array of knowns"; | ||
| 131 | |||
| 132 | case BackendDAE.MAIN(varData = BVariable.VAR_DATA_SIM(knowns = knowns)) | ||
| 133 | algorithm | ||
| 134 |
2/2✓ Branch 1 taken 1 time.
✓ Branch 2 taken 187 times.
|
188 | if Flags.isSet(Flags.JAC_DUMP) then |
| 135 | 1 | print(StringUtil.headline_1("[symjacdump] Creating symbolic Jacobians:") + "\n"); | |
| 136 | end if; | ||
| 137 | |||
| 138 | name := match kind | ||
| 139 | case NBPartition.Kind.ODE algorithm | ||
| 140 | name := "ODE_JAC"; | ||
| 141 | 187 | bdae.ode := applyToPartitions(bdae.ode, bdae.funcMap, knowns, name, func); | |
| 142 | then name; | ||
| 143 | case NBPartition.Kind.DAE algorithm | ||
| 144 | name := "DAE_JAC"; | ||
| 145 | 2 | bdae.dae := SOME(applyToPartitions(Util.getOption(bdae.dae), bdae.funcMap, knowns, name, func)); | |
| 146 | then name; | ||
| 147 | else algorithm | ||
| 148 | ✗ | Error.addMessage(Error.INTERNAL_ERROR,{getInstanceName() + " failed for: " + Partition.Partition.kindToString(kind)}); | |
| 149 | ✗ | then fail(); | |
| 150 | end match; | ||
| 151 | |||
| 152 | // DAE mode: SimCode reads only the DAE partition jacobian | ||
| 153 | 188 | bdae.ode_event := applyToPartitions(bdae.ode_event, bdae.funcMap, knowns, name, func, kind <> NBPartition.Kind.DAE); | |
| 154 | 188 | bdae.algebraic := applyToPartitions(bdae.algebraic, bdae.funcMap, knowns, name, func); | |
| 155 | 188 | bdae.alg_event := applyToPartitions(bdae.alg_event, bdae.funcMap, knowns, name, func); | |
| 156 | 188 | bdae.init := applyToPartitions(bdae.init, bdae.funcMap, knowns, name, func); | |
| 157 |
3/4✗ Branch 0 not taken.
✓ Branch 1 taken 188 times.
✓ Branch 2 taken 6 times.
✓ Branch 3 taken 182 times.
|
188 | if isSome(bdae.init_0) then |
| 158 | 12 | bdae.init_0 := SOME(applyToPartitions(Util.getOption(bdae.init_0), bdae.funcMap, knowns, name, func)); | |
| 159 | end if; | ||
| 160 | then bdae; | ||
| 161 | |||
| 162 | else algorithm | ||
| 163 | // maybe add failtrace here and allow failing | ||
| 164 | ✗ | Error.addMessage(Error.INTERNAL_ERROR,{getInstanceName() + " failed for: " + BackendDAE.toString(bdae)}); | |
| 165 | ✗ | then fail(); | |
| 166 | |||
| 167 | end match; | ||
| 168 | end main; | ||
| 169 | |||
| 170 | function applyToPartitions | ||
| 171 | input output list<Partition.Partition> partitions; | ||
| 172 | input output UnorderedMap<Path, Function> funcMap; | ||
| 173 | input VariablePointers knowns; | ||
| 174 | input String name; | ||
| 175 | input Module.jacobianInterface func; | ||
| 176 | input Boolean simJacobian = true "also create the partition jacobian"; | ||
| 177 | algorithm | ||
| 178 |
5/6✓ Branch 0 taken 489 times.
✓ Branch 1 taken 946 times.
✓ Branch 2 taken 489 times.
✓ Branch 3 taken 946 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 946 times.
|
1435 | partitions := list(partJacobian(part, funcMap, knowns, name, func, simJacobian) for part in partitions); |
| 179 | end applyToPartitions; | ||
| 180 | |||
| 181 | function nonlinear | ||
| 182 | input VariablePointers seedCandidates; | ||
| 183 | input VariablePointers partialCandidates; | ||
| 184 | input EquationPointers equations; | ||
| 185 | input array<StrongComponent> comps; | ||
| 186 | input Option<Adjacency.Matrix> full; | ||
| 187 | input UnorderedMap<Path, Function> funcMap; | ||
| 188 | input String name; | ||
| 189 | input Boolean staticAsContinuous; | ||
| 190 | output Option<Jacobian> jacobian; | ||
| 191 | protected | ||
| 192 | constant Module.jacobianInterface func = if Flags.isSet(Flags.NLS_ANALYTIC_JACOBIAN) | ||
| 193 | then jacobianSymbolic | ||
| 194 | else jacobianNumeric; | ||
| 195 | algorithm | ||
| 196 | try | ||
| 197 |
3/6✗ Branch 0 not taken.
✓ Branch 1 taken 85 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✓ Branch 5 taken 26 times.
✓ Branch 6 taken 59 times.
|
170 | jacobian := func( |
| 198 | name = name, | ||
| 199 | jacType = JacobianType.NLS, | ||
| 200 | seedCandidates = seedCandidates, | ||
| 201 | partialCandidates = partialCandidates, | ||
| 202 | equations = equations, | ||
| 203 | strongComponents = SOME(comps), | ||
| 204 | full = full, | ||
| 205 | funcMap = funcMap, | ||
| 206 | staticAsContinuous = staticAsContinuous | ||
| 207 | ); | ||
| 208 | else | ||
| 209 | // not everything can be differentiated symbolically, e.g. functions with function inputs | ||
| 210 | 1 | jacobian := jacobianNumeric( | |
| 211 | name = name, | ||
| 212 | jacType = JacobianType.NLS, | ||
| 213 | seedCandidates = seedCandidates, | ||
| 214 | partialCandidates = partialCandidates, | ||
| 215 | equations = equations, | ||
| 216 | strongComponents = SOME(comps), | ||
| 217 | full = full, | ||
| 218 | funcMap = funcMap, | ||
| 219 | staticAsContinuous = staticAsContinuous | ||
| 220 | ); | ||
| 221 | end try; | ||
| 222 | end nonlinear; | ||
| 223 | |||
| 224 | function combine | ||
| 225 | input list<BackendDAE> jacobians; | ||
| 226 | input String name; | ||
| 227 | output BackendDAE jacobian; | ||
| 228 | protected | ||
| 229 | JacobianType jacType = JacobianType.NLS; | ||
| 230 | list<Pointer<Variable>> variables = {}, unknowns = {}, auxiliaryVars = {}, aliasVars = {}; | ||
| 231 | list<Pointer<Variable>> diffVars = {}, dependencies = {}, resultVars = {}, tmpVars = {}, seedVars = {}; | ||
| 232 | list<StrongComponent> comps = {}; | ||
| 233 | list<Adjacency.Matrix> sparsity_patterns = {}; | ||
| 234 | VarData varData; | ||
| 235 | algorithm | ||
| 236 |
2/2✓ Branch 1 taken 80 times.
✓ Branch 2 taken 2 times.
|
82 | if List.hasOneElement(jacobians) then |
| 237 | 80 | jacobian := listHead(jacobians); | |
| 238 | jacobian := match jacobian case BackendDAE.JACOBIAN() algorithm | ||
| 239 | 80 | jacobian.name := name; | |
| 240 | then jacobian; | ||
| 241 | else algorithm | ||
| 242 | ✗ | Error.addMessage(Error.INTERNAL_ERROR,{getInstanceName() + " failed for\n" + BackendDAE.toString(jacobian)}); | |
| 243 | ✗ | then fail(); | |
| 244 | end match; | ||
| 245 | else | ||
| 246 |
2/2✓ Branch 0 taken 4 times.
✓ Branch 1 taken 2 times.
|
6 | for jac in jacobians loop |
| 247 | () := match jac | ||
| 248 | local | ||
| 249 | VarData tmpVarData; | ||
| 250 | |||
| 251 | case BackendDAE.JACOBIAN(varData = tmpVarData as VarData.VAR_DATA_JAC()) algorithm | ||
| 252 | 4 | jacType := jac.jacType; | |
| 253 | 4 | variables := listAppend(VariablePointers.toList(tmpVarData.variables), variables); | |
| 254 | 4 | unknowns := listAppend(VariablePointers.toList(tmpVarData.unknowns), unknowns); | |
| 255 | 4 | auxiliaryVars := listAppend(VariablePointers.toList(tmpVarData.auxiliaries), auxiliaryVars); | |
| 256 | 4 | aliasVars := listAppend(VariablePointers.toList(tmpVarData.aliasVars), aliasVars); | |
| 257 | 4 | diffVars := listAppend(VariablePointers.toList(tmpVarData.diffVars), diffVars); | |
| 258 | 4 | dependencies := listAppend(VariablePointers.toList(tmpVarData.dependencies), dependencies); | |
| 259 | 4 | resultVars := listAppend(VariablePointers.toList(tmpVarData.resultVars), resultVars); | |
| 260 | 4 | tmpVars := listAppend(VariablePointers.toList(tmpVarData.tmpVars), tmpVars); | |
| 261 | 4 | seedVars := listAppend(VariablePointers.toList(tmpVarData.seedVars), seedVars); | |
| 262 | 4 | comps := listAppend(arrayList(jac.comps), comps); | |
| 263 | 4 | sparsity_patterns := jac.sparsity :: sparsity_patterns; | |
| 264 | then (); | ||
| 265 | |||
| 266 | else algorithm | ||
| 267 | ✗ | Error.addMessage(Error.INTERNAL_ERROR,{getInstanceName() + " failed for\n" + BackendDAE.toString(jac)}); | |
| 268 | ✗ | then fail(); | |
| 269 | end match; | ||
| 270 | end for; | ||
| 271 | |||
| 272 | 2 | varData := VarData.VAR_DATA_JAC( | |
| 273 | variables = VariablePointers.fromList(variables), | ||
| 274 | unknowns = VariablePointers.fromList(unknowns), | ||
| 275 | auxiliaries = VariablePointers.fromList(auxiliaryVars), | ||
| 276 | aliasVars = VariablePointers.fromList(aliasVars), | ||
| 277 | diffVars = VariablePointers.fromList(diffVars), | ||
| 278 | dependencies = VariablePointers.fromList(dependencies), | ||
| 279 | resultVars = VariablePointers.fromList(resultVars), | ||
| 280 | tmpVars = VariablePointers.fromList(tmpVars), | ||
| 281 | seedVars = VariablePointers.fromList(seedVars) | ||
| 282 | ); | ||
| 283 | |||
| 284 |
3/4✓ Branch 0 taken 1 time.
✓ Branch 1 taken 1 time.
✓ Branch 3 taken 1 time.
✗ Branch 4 not taken.
|
3 | jacobian := BackendDAE.JACOBIAN( |
| 285 | name = name, | ||
| 286 | jacType = jacType, | ||
| 287 | varData = varData, | ||
| 288 | comps = listArray(comps), | ||
| 289 | sparsity = Adjacency.Matrix.combine(sparsity_patterns), | ||
| 290 | isAdjoint = name == "ADJ" // this is maybe bad (e.g. when name changes) | ||
| 291 | ); | ||
| 292 | end if; | ||
| 293 | end combine; | ||
| 294 | |||
| 295 | function getModule | ||
| 296 | "Returns the module function that was chosen by the user." | ||
| 297 | output Module.jacobianInterface func; | ||
| 298 | algorithm | ||
| 299 | func := match Flags.getConfigString(Flags.GENERATE_DYNAMIC_JACOBIAN) | ||
| 300 | case "symbolic" then jacobianSymbolic; | ||
| 301 | case "symbolicadjoint" then jacobianSymbolicAdjoint; | ||
| 302 | case "bidirectional" then jacobianSymbolic; | ||
| 303 | case "numeric" then jacobianNumeric; | ||
| 304 | case "none" then jacobianNone; | ||
| 305 | else algorithm | ||
| 306 | ✗ | Error.addMessage(Error.INTERNAL_ERROR,{getInstanceName() + " failed because of unknown jacobian type: " + Flags.getConfigString(Flags.GENERATE_DYNAMIC_JACOBIAN)}); | |
| 307 | ✗ | then fail(); | |
| 308 | end match; | ||
| 309 | end getModule; | ||
| 310 | |||
| 311 | function toString | ||
| 312 | input BackendDAE jacobian; | ||
| 313 | input output String str; | ||
| 314 | algorithm | ||
| 315 | 4 | str := BackendDAE.toString(jacobian, str); | |
| 316 | end toString; | ||
| 317 | |||
| 318 | function jacobianTypeString | ||
| 319 | input JacobianType jacType; | ||
| 320 | output String str; | ||
| 321 | algorithm | ||
| 322 | str := match jacType | ||
| 323 | case JacobianType.ODE then "[ODE]"; | ||
| 324 | case JacobianType.DAE then "[DAE]"; | ||
| 325 | case JacobianType.LS then "[LS-]"; | ||
| 326 | case JacobianType.NLS then "[NLS]"; | ||
| 327 | case JacobianType.OPT_LFG then "[OPT-LFG]"; | ||
| 328 | case JacobianType.OPT_MRF then "[OPT-MRF]"; | ||
| 329 | case JacobianType.OPT_R0 then "[OPT-R0]"; | ||
| 330 | else "[ERR]"; | ||
| 331 | end match; | ||
| 332 | end jacobianTypeString; | ||
| 333 | |||
| 334 | |||
| 335 | protected | ||
| 336 | // ToDo: all the DAEMode stuff is probably incorrect! | ||
| 337 | |||
| 338 | // TODO: refactor with map | ||
| 339 | function getOptimizableVars | ||
| 340 | input VariablePointers variables; | ||
| 341 | output list<Pointer<Variable>> optimizable_vars = {}; | ||
| 342 | algorithm | ||
| 343 | ✗ | for var_ptr in VariablePointers.toList(variables) loop | |
| 344 | ✗ | if BVariable.isOptimizable(var_ptr) then | |
| 345 | optimizable_vars := var_ptr :: optimizable_vars; | ||
| 346 | end if; | ||
| 347 | end for; | ||
| 348 | end getOptimizableVars; | ||
| 349 | |||
| 350 | function getSeedCandidatesDynamicOptimization | ||
| 351 | input Partition.Partition part; | ||
| 352 | input VariablePointers all_knowns; | ||
| 353 | input BVariable.checkVar filter; | ||
| 354 | output list<Pointer<Variable>> unknowns; | ||
| 355 | protected | ||
| 356 | list<Pointer<Variable>> derivative_vars, unknown_states; | ||
| 357 | algorithm | ||
| 358 | // we could absorb the filter into getOptimizableVars as its faster | ||
| 359 | ✗ | unknowns := getOptimizableVars(all_knowns); // all optimizable inputs + parameters | |
| 360 | ✗ | derivative_vars := list(var for var guard(BVariable.isStateDerivative(var)) in VariablePointers.toList(part.unknowns)); | |
| 361 | ✗ | unknown_states := list(Util.getOption(BVariable.getVarState(var)) for var in derivative_vars); // all states | |
| 362 | ✗ | unknowns := listAppend(unknown_states, unknowns); // all states, inputs and parameters (optimizable) | |
| 363 | ✗ | unknowns := List.filterOnTrue(unknowns, filter); | |
| 364 | // sort? | ||
| 365 | end getSeedCandidatesDynamicOptimization; | ||
| 366 | |||
| 367 | function getLfgPartialCandidates | ||
| 368 | input Partition.Partition part; | ||
| 369 | output list<Pointer<Variable>> partialCandidates; | ||
| 370 | protected | ||
| 371 | list<Pointer<Variable>> lagrange_vars = {}, derivative_vars = {}, path_vars = {}; | ||
| 372 | algorithm | ||
| 373 | ✗ | for var_ptr in VariablePointers.toList(part.unknowns) loop | |
| 374 | ✗ | if BVariable.isLagrange(var_ptr) then | |
| 375 | lagrange_vars := var_ptr :: lagrange_vars; | ||
| 376 | elseif BVariable.isStateDerivative(var_ptr) then | ||
| 377 | derivative_vars := var_ptr :: derivative_vars; | ||
| 378 | elseif BVariable.isPathConstraint(var_ptr) then | ||
| 379 | path_vars := var_ptr :: path_vars; | ||
| 380 | end if; | ||
| 381 | end for; | ||
| 382 | ✗ | partialCandidates := listReverse(listAppend(lagrange_vars, listAppend(derivative_vars, path_vars))); | |
| 383 | end getLfgPartialCandidates; | ||
| 384 | |||
| 385 | function getMrfPartialCandidates | ||
| 386 | input Partition.Partition part; | ||
| 387 | output list<Pointer<Variable>> partialCandidates; | ||
| 388 | protected | ||
| 389 | list<Pointer<Variable>> mayer_vars = {}, final_vars = {}; | ||
| 390 | algorithm | ||
| 391 | ✗ | for var_ptr in VariablePointers.toList(part.unknowns) loop | |
| 392 | ✗ | if BVariable.isMayer(var_ptr) then | |
| 393 | mayer_vars := var_ptr :: mayer_vars; | ||
| 394 | elseif BVariable.isFinalConstraint(var_ptr) then | ||
| 395 | final_vars := var_ptr :: final_vars; | ||
| 396 | end if; | ||
| 397 | end for; | ||
| 398 | ✗ | partialCandidates := listReverse(listAppend(mayer_vars, final_vars)); | |
| 399 | end getMrfPartialCandidates; | ||
| 400 | |||
| 401 | function getR0PartialCandidates | ||
| 402 | input Partition.Partition part; | ||
| 403 | output list<Pointer<Variable>> partialCandidates = {}; | ||
| 404 | algorithm | ||
| 405 | ✗ | for var_ptr in VariablePointers.toList(part.unknowns) loop | |
| 406 | ✗ | if BVariable.isInitialConstraint(var_ptr) then | |
| 407 | partialCandidates := var_ptr :: partialCandidates; | ||
| 408 | end if; | ||
| 409 | end for; | ||
| 410 | ✗ | partialCandidates := listReverse(partialCandidates); | |
| 411 | end getR0PartialCandidates; | ||
| 412 | |||
| 413 | // TODO: before this is ever called, we should check if the variable / annotation pairs are even valid: e.g. path constraint with final time or so! | ||
| 414 | // add a module for optimization? where we check the model, may do some transformations etc? | ||
| 415 | function partJacobianDynamicOptimization | ||
| 416 | input Partition.Partition part; | ||
| 417 | input VariablePointers all_knowns; | ||
| 418 | input String name; | ||
| 419 | input Module.jacobianInterface func; | ||
| 420 | input UnorderedMap<Path, Function> funcMap; | ||
| 421 | output Option<Jacobian> LFG_jacobian; | ||
| 422 | output Option<Jacobian> MRF_jacobian; | ||
| 423 | output Option<Jacobian> R0_jacobian; | ||
| 424 | protected | ||
| 425 | Boolean staticAsContinuous = true; | ||
| 426 | VariablePointers seedCandidates, partialCandidates; | ||
| 427 | algorithm | ||
| 428 | // Lfg Jacobian (Lagrange (L), ODE (f), Path Constraints (g)), append all unkowns of partition, as we might need their partials for inner | ||
| 429 | ✗ | partialCandidates := VariablePointers.fromList(listAppend(getLfgPartialCandidates(part), VariablePointers.toList(part.unknowns)), part.unknowns.scalarized); | |
| 430 | ✗ | seedCandidates := VariablePointers.fromList(getSeedCandidatesDynamicOptimization(part, all_knowns, BVariable.isLfgVariable), partialCandidates.scalarized); | |
| 431 | |||
| 432 | // TODO: add _OPT to name? | ||
| 433 | ✗ | LFG_jacobian := func(name, JacobianType.OPT_LFG, seedCandidates, partialCandidates, | |
| 434 | part.equations, part.strongComponents, part.adjacencyMatrix, funcMap, staticAsContinuous); | ||
| 435 | |||
| 436 | // Mrf Jacobian (Mayer (M), Final Constraints (rf)), append all unkowns of partition, as we might need their partials | ||
| 437 | ✗ | partialCandidates := VariablePointers.fromList(listAppend(getMrfPartialCandidates(part), VariablePointers.toList(part.unknowns)), part.unknowns.scalarized); | |
| 438 | ✗ | seedCandidates := VariablePointers.fromList(getSeedCandidatesDynamicOptimization(part, all_knowns, BVariable.isMrfVariable), partialCandidates.scalarized); | |
| 439 | |||
| 440 | // TODO: add _OPT to name? | ||
| 441 | ✗ | MRF_jacobian := func(name, JacobianType.OPT_MRF, seedCandidates, partialCandidates, | |
| 442 | part.equations, part.strongComponents, part.adjacencyMatrix, funcMap, staticAsContinuous); | ||
| 443 | |||
| 444 | // r0 Jacobian (Initial Constraints (r0)), append all unkowns of partition, as we might need their partials | ||
| 445 | ✗ | partialCandidates := VariablePointers.fromList(listAppend(getR0PartialCandidates(part), VariablePointers.toList(part.unknowns)), part.unknowns.scalarized); | |
| 446 | ✗ | seedCandidates := VariablePointers.fromList(getSeedCandidatesDynamicOptimization(part, all_knowns, BVariable.isR0Variable), partialCandidates.scalarized); | |
| 447 | |||
| 448 | // TODO: add _OPT to name? | ||
| 449 | ✗ | R0_jacobian := func(name, JacobianType.OPT_R0, seedCandidates, partialCandidates, | |
| 450 | part.equations, part.strongComponents, part.adjacencyMatrix, funcMap, staticAsContinuous); | ||
| 451 | end partJacobianDynamicOptimization; | ||
| 452 | |||
| 453 | function partJacobian | ||
| 454 | input output Partition.Partition part; | ||
| 455 | input UnorderedMap<Path, Function> funcMap; | ||
| 456 | input VariablePointers knowns; | ||
| 457 | input String name "Context name for jacobian"; | ||
| 458 | input Module.jacobianInterface func; | ||
| 459 | input Boolean simJacobian = true; | ||
| 460 | protected | ||
| 461 | JacobianType jacType; | ||
| 462 | VariablePointers unknowns; | ||
| 463 | list<Pointer<Variable>> derivative_vars, state_vars; | ||
| 464 | VariablePointers seedCandidates, partialCandidates; | ||
| 465 | Option<Jacobian> jacobian, LFG_jacobian = NONE(), MRF_jacobian = NONE(), R0_jacobian = NONE() "Resulting jacobians"; | ||
| 466 | Option<Jacobian> adjointJac; | ||
| 467 | Partition.Kind kind = Partition.Partition.getKind(part); | ||
| 468 | Boolean updated; | ||
| 469 | algorithm | ||
| 470 | // create algebraic loop jacobians | ||
| 471 | part.strongComponents := match part.strongComponents | ||
| 472 | local | ||
| 473 | array<StrongComponent> comps; | ||
| 474 | StrongComponent tmp; | ||
| 475 | case SOME(comps) algorithm | ||
| 476 |
2/2✓ Branch 0 taken 487 times.
✓ Branch 1 taken 2 times.
|
6479 | for i in 1:arrayLength(comps) loop |
| 477 | 5990 | (tmp, updated) := compJacobian(comps[i], part.adjacencyMatrix, funcMap, kind); | |
| 478 |
2/2✓ Branch 0 taken 85 times.
✓ Branch 1 taken 5905 times.
|
5990 | if updated then arrayUpdate(comps, i, tmp); end if; |
| 479 | end for; | ||
| 480 | then SOME(comps); | ||
| 481 | else part.strongComponents; | ||
| 482 | end match; | ||
| 483 | |||
| 484 | // create the simulation jacobian | ||
| 485 |
3/4✗ Branch 0 not taken.
✓ Branch 1 taken 489 times.
✓ Branch 3 taken 406 times.
✓ Branch 4 taken 83 times.
|
489 | if simJacobian and Partition.Partition.isODEorDAE(part) then |
| 486 | 83 | partialCandidates := part.unknowns; | |
| 487 |
2/2✓ Branch 1 taken 1 time.
✓ Branch 2 taken 82 times.
|
83 | unknowns := if Partition.Partition.getKind(part) == NBPartition.Kind.DAE then Util.getOption(part.daeUnknowns) else part.unknowns; |
| 488 |
2/2✓ Branch 1 taken 82 times.
✓ Branch 2 taken 1 time.
|
83 | jacType := if Partition.Partition.getKind(part) == NBPartition.Kind.DAE then JacobianType.DAE else JacobianType.ODE; |
| 489 | |||
| 490 |
6/6✓ Branch 2 taken 425 times.
✓ Branch 3 taken 118 times.
✓ Branch 4 taken 543 times.
✓ Branch 5 taken 83 times.
✓ Branch 6 taken 118 times.
✓ Branch 7 taken 83 times.
|
626 | derivative_vars := list(var for var guard(BVariable.isStateDerivative(var)) in VariablePointers.toList(unknowns)); |
| 491 |
4/4✓ Branch 0 taken 118 times.
✓ Branch 1 taken 83 times.
✓ Branch 2 taken 118 times.
✓ Branch 3 taken 83 times.
|
201 | state_vars := list(Util.getOption(BVariable.getVarState(var)) for var in derivative_vars); |
| 492 | 83 | seedCandidates := VariablePointers.fromList(state_vars, partialCandidates.scalarized); | |
| 493 | |||
| 494 |
2/6✗ Branch 0 not taken.
✓ Branch 1 taken 83 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 7 taken 83 times.
✗ Branch 8 not taken.
|
166 | jacobian := func(name, jacType, seedCandidates, partialCandidates, part.equations, part.strongComponents, part.adjacencyMatrix, funcMap, Partition.kindIsInitial(kind)); |
| 495 | |||
| 496 |
1/2✗ Branch 1 not taken.
✓ Branch 2 taken 83 times.
|
83 | if Flags.getConfigBool(Flags.MOO_DYNAMIC_OPTIMIZATION) then |
| 497 | /* Add Lfg + Mr Jacobians for MOO dynamic optimization */ | ||
| 498 | ✗ | (LFG_jacobian, MRF_jacobian, R0_jacobian) := partJacobianDynamicOptimization(part, knowns, name, func, funcMap); | |
| 499 | end if; | ||
| 500 | |||
| 501 |
6/10✓ Branch 1 taken 1 time.
✓ Branch 2 taken 82 times.
✓ Branch 5 taken 1 time.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
✓ Branch 8 taken 1 time.
✓ Branch 9 taken 1 time.
✗ Branch 10 not taken.
✓ Branch 13 taken 1 time.
✗ Branch 14 not taken.
|
83 | if Flags.getConfigString(Flags.GENERATE_DYNAMIC_JACOBIAN) == "bidirectional" and isSome(jacobian) and not BackendDAE.getIsAdjoint(Util.getOption(jacobian)) then |
| 502 | // Bidirectional: generate adjoint jacobian in addition to forward | ||
| 503 | 1 | adjointJac := jacobianSymbolicAdjoint(name, jacType, seedCandidates, partialCandidates, part.equations, part.strongComponents, part.adjacencyMatrix, funcMap, kind == NBPartition.Kind.INI); | |
| 504 | 2 | part.association := Partition.Association.CONTINUOUS(kind, jacobian, adjointJac, LFG_jacobian, MRF_jacobian, R0_jacobian); | |
| 505 | elseif isSome(jacobian) then | ||
| 506 |
2/2✓ Branch 2 taken 3 times.
✓ Branch 3 taken 79 times.
|
82 | if BackendDAE.getIsAdjoint(Util.getOption(jacobian)) then |
| 507 | 6 | part.association := Partition.Association.CONTINUOUS(kind, NONE(), jacobian, LFG_jacobian, MRF_jacobian, R0_jacobian); | |
| 508 | else | ||
| 509 | 158 | part.association := Partition.Association.CONTINUOUS(kind, jacobian, NONE(), LFG_jacobian, MRF_jacobian, R0_jacobian); | |
| 510 | end if; | ||
| 511 | else | ||
| 512 | ✗ | part.association := Partition.Association.CONTINUOUS(kind, NONE(), NONE(), LFG_jacobian, MRF_jacobian, R0_jacobian); | |
| 513 | end if; | ||
| 514 |
2/2✓ Branch 1 taken 82 times.
✓ Branch 2 taken 1 time.
|
83 | if Flags.isSet(Flags.JAC_DUMP) then |
| 515 | 1 | print(Partition.Partition.toString(part, 2)); | |
| 516 | end if; | ||
| 517 | end if; | ||
| 518 | end partJacobian; | ||
| 519 | |||
| 520 | function forEquationStart | ||
| 521 | "Returns the INTEGER start value of a FOR_EQUATION's single iterator range. | ||
| 522 | Returns 0 when the equation is not a FOR_EQUATION or start is not an integer." | ||
| 523 | input Equation eqn; | ||
| 524 | output Integer s = 0; | ||
| 525 | protected | ||
| 526 | Expression start_exp; | ||
| 527 | algorithm | ||
| 528 | () := match eqn | ||
| 529 | case Equation.FOR_EQUATION(iter = Iterator.SINGLE(range = Expression.RANGE(start = start_exp))) | ||
| 530 | 9 | algorithm s := Expression.integerValueOrDefault(start_exp, 0); then (); | |
| 531 | else (); | ||
| 532 | end match; | ||
| 533 | end forEquationStart; | ||
| 534 | |||
| 535 | function withInnerComps | ||
| 536 | "an algebraic loop preceded by its inner components" | ||
| 537 | input StrongComponent comp; | ||
| 538 | output list<StrongComponent> comps; | ||
| 539 | algorithm | ||
| 540 | comps := match comp | ||
| 541 | local | ||
| 542 | Tearing strict; | ||
| 543 | 2 | case StrongComponent.ALGEBRAIC_LOOP(strict = strict) then listAppend(arrayList(strict.innerEquations), {comp}); | |
| 544 | else {comp}; | ||
| 545 | end match; | ||
| 546 | end withInnerComps; | ||
| 547 | |||
| 548 | function partialSliceSeedCandidates | ||
| 549 | "Creates per-element seed candidates for partial iteration-var slices where | ||
| 550 | the for-loop starts at or above the slice's first element index. This | ||
| 551 | avoids phantom seeds (e.g. $SEED.x[1] when x[2..N] is the NLS slice) that | ||
| 552 | would produce a zero sparsity column and cause sparsitySanityCheck to fail. | ||
| 553 | Falls back to whole-array seeds when the for-loop starts below the slice | ||
| 554 | (e.g. slice_for where i1 starts at 1 but x is sliced from x[2])." | ||
| 555 | input list<Slice<VariablePointer>> iteration_vars; | ||
| 556 | input list<Slice<BEquation.EquationPointer>> residual_eqns; | ||
| 557 | output list<VariablePointer> seed_candidates = {}; | ||
| 558 | protected | ||
| 559 | Integer for_start = 0; // 0 = no FOR_EQUATION found | ||
| 560 | Integer s; | ||
| 561 | Integer slice_first_1based; | ||
| 562 | list<Variable> elem_vars; | ||
| 563 | Variable var_elem; | ||
| 564 | algorithm | ||
| 565 | // Find minimum for-loop start across all FOR_EQUATION residuals. | ||
| 566 |
2/2✓ Branch 0 taken 727 times.
✓ Branch 1 taken 85 times.
|
812 | for eqn_slice in residual_eqns loop |
| 567 | 727 | s := forEquationStart(Pointer.access(Slice.getT(eqn_slice))); | |
| 568 |
2/2✓ Branch 0 taken 9 times.
✓ Branch 1 taken 718 times.
|
727 | if s > 0 then |
| 569 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 9 times.
|
9 | if for_start == 0 then |
| 570 | for_start := s; | ||
| 571 | else | ||
| 572 | for_start := intMin(for_start, s); | ||
| 573 | end if; | ||
| 574 | end if; | ||
| 575 | end for; | ||
| 576 | // Build seed candidates: per-element for safe partial slices, whole-array otherwise. | ||
| 577 |
2/2✓ Branch 0 taken 645 times.
✓ Branch 1 taken 85 times.
|
730 | for var_slice in iteration_vars loop |
| 578 |
2/2✓ Branch 0 taken 611 times.
✓ Branch 1 taken 34 times.
|
645 | if listEmpty(var_slice.indices) then |
| 579 | // Full slice: no phantom risk. | ||
| 580 | 611 | seed_candidates := Slice.getT(var_slice) :: seed_candidates; | |
| 581 | else | ||
| 582 | // Partial slice: first 0-based index + 1 gives the 1-based start element. | ||
| 583 | 34 | var_elem := Pointer.access(Slice.getT(var_slice)); | |
| 584 | 34 | slice_first_1based := listHead(var_slice.indices) + 1; | |
| 585 |
2/6✓ Branch 0 taken 34 times.
✗ Branch 1 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 34 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
|
34 | if (for_start == 0 or for_start >= slice_first_1based) and |
| 586 | (Type.isReal(Type.arrayElementType(var_elem.ty)) or Type.isComplex(Type.arrayElementType(var_elem.ty))) then | ||
| 587 | // Safe to use per-element seeds: the loop never reaches below the slice start | ||
| 588 | // (for_start >= slice_first_1based, or no FOR_EQUATION residual at all). | ||
| 589 | // Record (Complex) element types are included here too: scalarizeBackendVariable | ||
| 590 | // only splits the ARRAY dimension (still record-typed per-element results), and | ||
| 591 | // the later general VariablePointers.scalarize pass (already used for every other | ||
| 592 | // Jacobian's seedVars, e.g. NSimJacobian.mo) recurses into scalarizeComplexVariable | ||
| 593 | // to flatten those into scalar leaf fields - the same path record-typed FULL-slice | ||
| 594 | // seed candidates already go through above. Restricting to Real here (as before) | ||
| 595 | // meant a record-element partial slice fell back to the WHOLE parent array/record | ||
| 596 | // as a single seed candidate, which the general scalarize pass then expands to far | ||
| 597 | // more scalar columns than this Jacobian's true torn-unknown count (sizeCols ended | ||
| 598 | // up way bigger than the NLS's own size), silently degrading every such Jacobian to | ||
| 599 | // a numeric one at runtime (see nonlinearSystem.c's sizeCols-vs-size safety check). | ||
| 600 | // NOTE: for_start was computed above but previously never consulted here (only | ||
| 601 | // "slice_first_1based > 1" was checked) -- a dead-code bug that silently defeated | ||
| 602 | // exactly the case this function's own docstring describes (slice_for.mos: the | ||
| 603 | // for-loop starts at 1 but x is sliced from x[2], so a symbolic body term like | ||
| 604 | // x[$i1] at $i1=1 needs a seed for x[1], which per-element scalarization of just | ||
| 605 | // x[2..4] can never provide). | ||
| 606 | // the slice indices refer to the resized sizes of resizable dimensions | ||
| 607 | 34 | elem_vars := Scalarize.scalarizeBackendVariable(var_elem, var_slice.indices, resize = true); | |
| 608 |
2/2✓ Branch 0 taken 122 times.
✓ Branch 1 taken 34 times.
|
156 | for v in elem_vars loop |
| 609 | 122 | seed_candidates := Pointer.create(v) :: seed_candidates; | |
| 610 | end for; | ||
| 611 | else | ||
| 612 | // Unsafe or no phantom risk avoidance possible: fall back to whole-array seed pointer. | ||
| 613 | ✗ | seed_candidates := Slice.getT(var_slice) :: seed_candidates; | |
| 614 | end if; | ||
| 615 | end if; | ||
| 616 | end for; | ||
| 617 | 85 | seed_candidates := listReverse(seed_candidates); | |
| 618 | end partialSliceSeedCandidates; | ||
| 619 | |||
| 620 | function compJacobian | ||
| 621 | input output StrongComponent comp; | ||
| 622 | input Option<Adjacency.Matrix> full; | ||
| 623 | input UnorderedMap<Path, Function> funcMap; | ||
| 624 | input Partition.Kind kind; | ||
| 625 | output Boolean updated; | ||
| 626 | protected | ||
| 627 | Tearing strict; | ||
| 628 | list<StrongComponent> residual_comps; | ||
| 629 | list<VariablePointer> seed_candidates, residual_vars, inner_vars; | ||
| 630 | constant Boolean staticAsContinuous = Partition.kindIsInitial(kind); | ||
| 631 | algorithm | ||
| 632 | (comp, updated) := match comp | ||
| 633 | // nothing to differentiate if all iteration variables are discrete (e.g. Boolean) | ||
| 634 | case StrongComponent.ALGEBRAIC_LOOP(strict = strict) | ||
| 635 | guard(not List.any(list(Slice.getT(v) for v in strict.iteration_vars), function BVariable.isContinuous(staticAsContinuous = staticAsContinuous))) | ||
| 636 | then (comp, false); | ||
| 637 | |||
| 638 | case StrongComponent.ALGEBRAIC_LOOP(strict = strict) algorithm | ||
| 639 | // create residual components | ||
| 640 |
4/4✓ Branch 0 taken 727 times.
✓ Branch 1 taken 85 times.
✓ Branch 2 taken 727 times.
✓ Branch 3 taken 85 times.
|
812 | residual_comps := list(StrongComponent.fromSolvedEquationSlice(eqn) for eqn in strict.residual_eqns); |
| 641 | |||
| 642 | // create seed and partial candidates | ||
| 643 | 85 | seed_candidates := partialSliceSeedCandidates(strict.iteration_vars, strict.residual_eqns); | |
| 644 |
4/4✓ Branch 0 taken 727 times.
✓ Branch 1 taken 85 times.
✓ Branch 2 taken 727 times.
✓ Branch 3 taken 85 times.
|
812 | residual_vars := list(Equation.getResidualVar(Slice.getT(eqn)) for eqn in strict.residual_eqns); |
| 645 |
10/10✓ Branch 0 taken 403 times.
✓ Branch 1 taken 85 times.
✓ Branch 3 taken 403 times.
✓ Branch 4 taken 85 times.
✓ Branch 8 taken 5 times.
✓ Branch 9 taken 401 times.
✓ Branch 10 taken 406 times.
✓ Branch 11 taken 403 times.
✓ Branch 12 taken 401 times.
✓ Branch 13 taken 403 times.
|
979 | inner_vars := listAppend(list(var for var guard(BVariable.isContinuous(var, staticAsContinuous)) in StrongComponent.getVariables(comp)) for comp in strict.innerEquations); |
| 646 | |||
| 647 | // update jacobian to take slices (just to have correct inner variables and such) | ||
| 648 |
8/8✓ Branch 0 taken 727 times.
✓ Branch 1 taken 85 times.
✓ Branch 2 taken 727 times.
✓ Branch 3 taken 85 times.
✓ Branch 5 taken 75 times.
✓ Branch 6 taken 10 times.
✓ Branch 7 taken 46 times.
✓ Branch 8 taken 29 times.
|
943 | strict.jac := nonlinear( |
| 649 | seedCandidates = VariablePointers.fromList(seed_candidates), | ||
| 650 | partialCandidates = VariablePointers.fromList(listAppend(residual_vars, inner_vars)), | ||
| 651 | equations = EquationPointers.fromList(list(Slice.getT(eqn) for eqn in strict.residual_eqns)), | ||
| 652 | comps = Array.appendList(strict.innerEquations, residual_comps), | ||
| 653 | full = full, | ||
| 654 | funcMap = funcMap, | ||
| 655 | // comp.idx is assigned from two disjoint, independently-numbered sources: | ||
| 656 | // NBSolve.mo's Tearing.implicit() (single-equation loops promoted from an | ||
| 657 | // unsolvable explicit equation) and NBTearing.mo's initialize() (genuine | ||
| 658 | // multi-equation torn loops), each restarting its own count at 1 per | ||
| 659 | // partition kind. Two components from different sources can therefore | ||
| 660 | // share the same (kind, idx), which previously collided on the exact same | ||
| 661 | // generated Jacobian name (e.g. both "INI_NLS_JAC_1") -- note comp.status | ||
| 662 | // alone can't distinguish them, since NBSolve.mo's solveStrongComponent | ||
| 663 | // marks EVERY algebraic loop status=IMPLICIT, torn or not. Tag | ||
| 664 | // comp.implicitlyCreated components distinctly so the two numbering | ||
| 665 | // spaces can never collide, without changing either one's actual numbers | ||
| 666 | // (which existing reference test outputs depend on). | ||
| 667 | name = Partition.Partition.kindToString(kind) | ||
| 668 | + (if comp.implicitlyCreated then "_NLS_IMPL_" | ||
| 669 | else if comp.linear then "_LS_JAC_" else "_NLS_JAC_") | ||
| 670 | + intString(comp.idx), | ||
| 671 | staticAsContinuous = staticAsContinuous); | ||
| 672 | 85 | comp.strict := strict; | |
| 673 | |||
| 674 |
2/2✓ Branch 1 taken 2 times.
✓ Branch 2 taken 83 times.
|
85 | if Flags.isSet(Flags.JAC_DUMP) then |
| 675 | 2 | print(StrongComponent.toString(comp) + "\n"); | |
| 676 | end if; | ||
| 677 | then (comp, true); | ||
| 678 | else (comp, false); | ||
| 679 | end match; | ||
| 680 | end compJacobian; | ||
| 681 | |||
| 682 | function jacobianSymbolic extends Module.jacobianInterface; | ||
| 683 | protected | ||
| 684 | list<StrongComponent> comps, diffed_comps; | ||
| 685 | Pointer<list<Pointer<Variable>>> seed_vars_ptr = Pointer.create({}); | ||
| 686 | Pointer<list<Pointer<Variable>>> pDer_vars_ptr = Pointer.create({}); | ||
| 687 | UnorderedMap<ComponentRef,ComponentRef> diff_map = UnorderedMap.new<ComponentRef>(ComponentRef.hash, ComponentRef.isEqual); | ||
| 688 | UnorderedMap<ComponentRef,ComponentRef> seed_diff_map; | ||
| 689 | Differentiate.DifferentiationArguments diffArguments; | ||
| 690 | Pointer<Integer> idx = Pointer.create(0); | ||
| 691 | |||
| 692 | VariablePointers adjacencyVars; | ||
| 693 | list<StrongComponent> sparsity_comps; | ||
| 694 | list<Pointer<Variable>> all_vars, unknown_vars, aux_vars, alias_vars, depend_vars, res_vars, res_vars_d, tmp_vars, tmp_vars_d, seed_vars, seed_vars_d; | ||
| 695 | BVariable.VarData varDataJac; | ||
| 696 | Adjacency.Matrix fullLocal, sparsity; | ||
| 697 | UnorderedSet<ComponentRef> seed_set = UnorderedSet.new(ComponentRef.hash, ComponentRef.isEqual); | ||
| 698 | UnorderedSet<ComponentRef> pder_set = UnorderedSet.new(ComponentRef.hash, ComponentRef.isEqual); | ||
| 699 | UnorderedSet<ComponentRef> adj_base_seen; | ||
| 700 | list<Pointer<Variable>> adj_seed_list; | ||
| 701 | ComponentRef adj_base_cref; | ||
| 702 | BVariable.checkVar func = getTmpFilterFunction(jacType); | ||
| 703 | algorithm | ||
| 704 |
2/4✗ Branch 0 not taken.
✓ Branch 1 taken 86 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 86 times.
|
86 | if isSome(strongComponents) then |
| 705 | // filter all discrete strong components and differentiate the others | ||
| 706 | // todo: mixed algebraic loops should be here without the discrete subsets | ||
| 707 |
6/6✓ Branch 3 taken 7 times.
✓ Branch 4 taken 1131 times.
✓ Branch 5 taken 1138 times.
✓ Branch 6 taken 86 times.
✓ Branch 7 taken 1131 times.
✓ Branch 8 taken 86 times.
|
2441 | comps := list(comp for comp guard(not StrongComponent.isDiscrete(comp)) in Util.getOption(strongComponents)); |
| 708 | else | ||
| 709 | ✗ | Error.addMessage(Error.INTERNAL_ERROR,{getInstanceName() + " failed because no strong components were given!"}); | |
| 710 | ✗ | fail(); | |
| 711 | end if; | ||
| 712 | |||
| 713 | // create seed vars | ||
| 714 |
2/2✓ Branch 0 taken 28 times.
✓ Branch 1 taken 58 times.
|
114 | VariablePointers.mapPtr(seedCandidates, function makeVarTraverse(name = name, vars_ptr = seed_vars_ptr, map = diff_map, |
| 715 | makeVar = BVariable.makeSeedVar, staticAsContinuous = staticAsContinuous)); | ||
| 716 |
2/2✓ Branch 1 taken 738 times.
✓ Branch 2 taken 86 times.
|
824 | for v in VariablePointers.toList(seedCandidates) loop |
| 717 |
1/2✓ Branch 1 taken 738 times.
✗ Branch 2 not taken.
|
738 | if BVariable.isContinuous(v, staticAsContinuous) then |
| 718 | 738 | UnorderedSet.add(BVariable.getVarName(v), seed_set); | |
| 719 | // Also add base cref so iterator-subscripted deps from for-loop equations | ||
| 720 | // (where subscript is an iterator variable, not a literal integer) can | ||
| 721 | // match via base fallback in filterSet. | ||
| 722 | 738 | UnorderedSet.add(ComponentRef.stripSubscriptsAll(BVariable.getVarName(v)), seed_set); | |
| 723 | end if; | ||
| 724 | end for; | ||
| 725 | |||
| 726 | // create pDer vars (also filters out discrete vars) | ||
| 727 | 86 | (res_vars, tmp_vars) := List.splitOnTrue(VariablePointers.toList(partialCandidates), func); | |
| 728 | 86 | (tmp_vars, _) := List.splitOnTrue(tmp_vars, function BVariable.isContinuous(staticAsContinuous = staticAsContinuous)); | |
| 729 | |||
| 730 |
2/2✓ Branch 0 taken 732 times.
✓ Branch 1 taken 86 times.
|
818 | for v in res_vars loop |
| 731 | 732 | makeVarTraverse(v, name, pDer_vars_ptr, diff_map, function BVariable.makePDerVar(isTmp = false), staticAsContinuous = staticAsContinuous); | |
| 732 | end for; | ||
| 733 | |||
| 734 |
2/2✓ Branch 0 taken 732 times.
✓ Branch 1 taken 86 times.
|
818 | for v in res_vars loop |
| 735 | 732 | UnorderedSet.add(BVariable.getVarName(v), pder_set); | |
| 736 | end for; | ||
| 737 | 86 | res_vars_d := listReverse(Pointer.access(pDer_vars_ptr)); | |
| 738 | |||
| 739 | 86 | pDer_vars_ptr := Pointer.create({}); | |
| 740 | // Snapshot diff_map before adding inner LS tmp pder entries. | ||
| 741 | // When an outer iter var and an inner LS var share the same base ComponentRef | ||
| 742 | // (slices of the same array variable), the tmp pder pass below would overwrite | ||
| 743 | // the outer seed entry. fullToSparsity must see the pre-overwrite version so | ||
| 744 | // that outer iter var dependencies resolve to outer seed columns (0..N-1), not | ||
| 745 | // to inner LS tmp pder columns (N..N+M-1). | ||
| 746 | 86 | seed_diff_map := UnorderedMap.copy(diff_map); | |
| 747 |
2/2✓ Branch 3 taken 311 times.
✓ Branch 4 taken 86 times.
|
397 | for v in tmp_vars loop makeVarTraverse(v, name, pDer_vars_ptr, diff_map, function BVariable.makePDerVar(isTmp = true), staticAsContinuous = staticAsContinuous); end for; |
| 748 | 86 | tmp_vars_d := Pointer.access(pDer_vars_ptr); | |
| 749 | |||
| 750 | // Build differentiation argument structure | ||
| 751 | 86 | diffArguments := Differentiate.DIFFERENTIATION_ARGUMENTS( | |
| 752 | diffCref = ComponentRef.EMPTY(), // no explicit cref necessary, rules are set by diff map | ||
| 753 | new_vars = {}, | ||
| 754 | diff_map = SOME(diff_map), // seed and temporary cref map | ||
| 755 | diffType = NBDifferentiate.DifferentiationType.JACOBIAN, | ||
| 756 | funcMap = funcMap, | ||
| 757 | scalarized = seedCandidates.scalarized, | ||
| 758 | adjoint_map = NONE(), | ||
| 759 | current_grad = Expression.EMPTY(Type.REAL()), | ||
| 760 | collectAdjoints = false | ||
| 761 | ); | ||
| 762 | |||
| 763 | // differentiate all strong components | ||
| 764 | 86 | (diffed_comps, diffArguments) := Differentiate.differentiateStrongComponentList(comps, diffArguments, idx, name, getInstanceName()); | |
| 765 | |||
| 766 | // collect var data (most of this can be removed) | ||
| 767 | 85 | unknown_vars := listAppend(res_vars_d, tmp_vars_d); | |
| 768 | all_vars := unknown_vars; // add other vars later on | ||
| 769 | |||
| 770 | 85 | seed_vars_d := listReverse(Pointer.access(seed_vars_ptr)); | |
| 771 | aux_vars := seed_vars_d; // add other auxiliaries later on | ||
| 772 | alias_vars := {}; | ||
| 773 | depend_vars := {}; | ||
| 774 | |||
| 775 | 85 | varDataJac := BVariable.VAR_DATA_JAC( | |
| 776 | variables = VariablePointers.fromList(all_vars), | ||
| 777 | unknowns = VariablePointers.fromList(unknown_vars), | ||
| 778 | auxiliaries = VariablePointers.fromList(aux_vars), | ||
| 779 | aliasVars = VariablePointers.fromList(alias_vars), | ||
| 780 | diffVars = partialCandidates, | ||
| 781 | dependencies = VariablePointers.fromList(depend_vars), | ||
| 782 | resultVars = VariablePointers.fromList(res_vars_d), | ||
| 783 | tmpVars = VariablePointers.fromList(tmp_vars_d), | ||
| 784 | seedVars = VariablePointers.fromList(seed_vars_d) | ||
| 785 | ); | ||
| 786 | |||
| 787 | // Always rebuild the full matrix from the actual comps equations. | ||
| 788 | // Using part.adjacencyMatrix (the `full` param) would fail because residual | ||
| 789 | // equations created by finalize() during tearing get new names and don't | ||
| 790 | // appear in the partition's pre-tearing adjacency matrix. | ||
| 791 | // Build adjacencyVars from unique base variable ptrs derived from seedCandidates. | ||
| 792 | // When seedCandidates contains scalar element ptrs for partial-slice NLS iter vars, | ||
| 793 | // all elements of the same array share the same base ptr. Using base ptrs here | ||
| 794 | // preserves pseudo=true subscript-stripped lookup in getDependentCref, which matches | ||
| 795 | // any element expression (e.g. module[i].T for iterator i) to the base column. | ||
| 796 | 85 | adj_base_seen := UnorderedSet.new(ComponentRef.hash, ComponentRef.isEqual); | |
| 797 | adj_seed_list := {}; | ||
| 798 |
2/2✓ Branch 1 taken 737 times.
✓ Branch 2 taken 85 times.
|
822 | for v in VariablePointers.toList(seedCandidates) loop |
| 799 | 737 | adj_base_cref := ComponentRef.stripSubscriptsAll(BVariable.getVarName(v)); | |
| 800 |
2/2✓ Branch 1 taken 649 times.
✓ Branch 2 taken 88 times.
|
737 | if not UnorderedSet.contains(adj_base_cref, adj_base_seen) then |
| 801 | 649 | UnorderedSet.add(adj_base_cref, adj_base_seen); | |
| 802 | 649 | adj_seed_list := BVariable.getVarPointer(BVariable.getVarName(v), sourceInfo()) :: adj_seed_list; | |
| 803 | end if; | ||
| 804 | end for; | ||
| 805 | 85 | adjacencyVars := VariablePointers.fromList(listReverse(adj_seed_list)); | |
| 806 | 85 | adjacencyVars := VariablePointers.addList(tmp_vars, adjacencyVars); | |
| 807 | // For ODE Jacobians, also include state derivatives as adjacency variables. | ||
| 808 | // Some equations use der(x_j) as an RHS input (e.g. der(x_i) = f(der(x_j), x_k)). | ||
| 809 | // Without this, the transitive seed dependency der(x_i) -> der(x_j) -> x_j is lost. | ||
| 810 |
2/2✓ Branch 0 taken 3 times.
✓ Branch 1 taken 82 times.
|
85 | if jacType == JacobianType.ODE then |
| 811 | 3 | adjacencyVars := VariablePointers.addList(res_vars, adjacencyVars); | |
| 812 | end if; | ||
| 813 | // with resizable arrays the inner equations of algebraic loops are needed for their dependencies | ||
| 814 |
6/6✓ Branch 1 taken 12 times.
✓ Branch 2 taken 73 times.
✓ Branch 3 taken 76 times.
✓ Branch 4 taken 12 times.
✓ Branch 5 taken 76 times.
✓ Branch 6 taken 12 times.
|
161 | sparsity_comps := if Flags.getConfigBool(Flags.RESIZABLE_ARRAYS) then List.flatten(list(withInnerComps(comp) for comp in comps)) else comps; |
| 815 |
4/4✓ Branch 0 taken 1130 times.
✓ Branch 1 taken 85 times.
✓ Branch 2 taken 1130 times.
✓ Branch 3 taken 85 times.
|
1215 | fullLocal := Adjacency.Matrix.createFull(adjacencyVars, |
| 816 | EquationPointers.fromList(List.flatten(list(StrongComponent.getEquations(comp) for comp in sparsity_comps)))); | ||
| 817 | 85 | sparsity := Adjacency.Matrix.fullToSparsity(fullLocal, sparsity_comps, seed_set, pder_set, seed_diff_map); | |
| 818 | |||
| 819 | 85 | jacobian := SOME(Jacobian.JACOBIAN( | |
| 820 | name = name, | ||
| 821 | jacType = jacType, | ||
| 822 | varData = varDataJac, | ||
| 823 | comps = listArray(diffed_comps), | ||
| 824 | sparsity = sparsity, | ||
| 825 | isAdjoint = false | ||
| 826 | )); | ||
| 827 | end jacobianSymbolic; | ||
| 828 | |||
| 829 | function sizeClassificationFromType | ||
| 830 | input Type ty; | ||
| 831 | output SizeClassification sc; | ||
| 832 | algorithm | ||
| 833 | sc := match Type.dimensionCount(ty) | ||
| 834 | case 0 then SizeClassification.SCALAR; | ||
| 835 | case 1 then SizeClassification.ELEMENT_WISE; | ||
| 836 | case 2 then SizeClassification.MATRIX; | ||
| 837 | else SizeClassification.ELEMENT_WISE; | ||
| 838 | end match; | ||
| 839 | end sizeClassificationFromType; | ||
| 840 | |||
| 841 | // Helper: build addition (or single term) expression from a list of terms for a given LHS cref. | ||
| 842 | function buildAdjointRhs | ||
| 843 | input ComponentRef lhsCref; | ||
| 844 | input list<Expression> terms; | ||
| 845 | output Expression rhs; | ||
| 846 | protected | ||
| 847 | Type vty; | ||
| 848 | SizeClassification sc; | ||
| 849 | Operator addOp; | ||
| 850 | algorithm | ||
| 851 | // Retrieve variable type | ||
| 852 | 3 | vty := ComponentRef.getComponentType(lhsCref); | |
| 853 | |||
| 854 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 3 times.
|
3 | if listEmpty(terms) then |
| 855 | ✗ | rhs := Expression.makeZero(vty); | |
| 856 | ✗ | return; | |
| 857 | end if; | ||
| 858 | |||
| 859 |
1/2✓ Branch 1 taken 3 times.
✗ Branch 2 not taken.
|
3 | if List.hasOneElement(terms) then |
| 860 | 3 | rhs := listHead(terms); | |
| 861 | 3 | return; | |
| 862 | end if; | ||
| 863 | |||
| 864 | ✗ | sc := sizeClassificationFromType(vty); | |
| 865 | ✗ | addOp := Operator.fromClassification( | |
| 866 | (MathClassification.ADDITION, sc), | ||
| 867 | vty | ||
| 868 | ); | ||
| 869 | |||
| 870 | ✗ | rhs := SimplifyExp.simplify(Expression.MULTARY(terms, {}, addOp)); | |
| 871 | ✗ | rhs := Expression.map(rhs, Expression.repairOperator); | |
| 872 | end buildAdjointRhs; | ||
| 873 | |||
| 874 | // Helper: run reverse-mode on a residual expression with a given seed (current_grad), | ||
| 875 | // accumulating into the provided adjoint_map. Returns updated DifferentiationArguments. | ||
| 876 | function accumulateAdjointForResidual | ||
| 877 | input Expression residual; | ||
| 878 | input Expression seed; // current_grad, typically a lambda_i cref | ||
| 879 | input UnorderedMap<ComponentRef,ComponentRef> diff_map; | ||
| 880 | input UnorderedMap<Path, Function> funcMapIn; | ||
| 881 | input Boolean scalarized; | ||
| 882 | input UnorderedMap<ComponentRef, AdjointTermList> adjoint_map_in; | ||
| 883 | output Differentiate.DifferentiationArguments diffArguments; | ||
| 884 | algorithm | ||
| 885 | // Prepare args to collect adjoints into the incoming map | ||
| 886 |
1/2✓ Branch 0 taken 1 time.
✗ Branch 1 not taken.
|
2 | diffArguments := Differentiate.DIFFERENTIATION_ARGUMENTS( |
| 887 | diffCref = ComponentRef.EMPTY(), | ||
| 888 | new_vars = {}, | ||
| 889 | diff_map = SOME(diff_map), | ||
| 890 | diffType = NBDifferentiate.DifferentiationType.JACOBIAN, | ||
| 891 | funcMap = funcMapIn, | ||
| 892 | scalarized = scalarized, | ||
| 893 | adjoint_map = SOME(adjoint_map_in), | ||
| 894 | current_grad = seed, | ||
| 895 | collectAdjoints = true | ||
| 896 | ); | ||
| 897 | |||
| 898 | // Run reverse-mode on the residual expression. | ||
| 899 | 1 | (_, diffArguments) := NBDifferentiate.differentiateExpression(residual, diffArguments); | |
| 900 | end accumulateAdjointForResidual; | ||
| 901 | |||
| 902 | // Reusable builder for a SINGLE_COMPONENT adjoint assignment (tmp or result var). | ||
| 903 | function makeAdjointComponentFromRhs | ||
| 904 | input ComponentRef lhsKey; | ||
| 905 | input Expression rhsExpr; | ||
| 906 | input String contextName; | ||
| 907 | input Integer eqIndex; | ||
| 908 | output NBStrongComponent diffed_comp; | ||
| 909 | protected | ||
| 910 | Pointer<NBEquation.Equation> eqPtr; | ||
| 911 | NBEquation.Equation eq; | ||
| 912 | Pointer<Variable> lhsVarPtr; | ||
| 913 | algorithm | ||
| 914 | ✗ | eqPtr := Equation.makeAssignment( | |
| 915 | Expression.fromCref(lhsKey), | ||
| 916 | rhsExpr, | ||
| 917 | Pointer.create(eqIndex), | ||
| 918 | contextName, | ||
| 919 | BEquation.Iterator.EMPTY(), | ||
| 920 | NBEquation.EquationAttributes.default(NBEquation.EquationKind.CONTINUOUS, false) | ||
| 921 | ); | ||
| 922 | |||
| 923 | ✗ | lhsVarPtr := BVariable.getVarPointer(lhsKey, sourceInfo()); | |
| 924 | ✗ | eq := Pointer.access(eqPtr); | |
| 925 | |||
| 926 | diffed_comp := match eq | ||
| 927 | case NBEquation.SCALAR_EQUATION() algorithm | ||
| 928 | ✗ | if not listEmpty(ComponentRef.subscriptsAllFlat(lhsKey)) then | |
| 929 | // Represent as a sliced component of size 1 | ||
| 930 | ✗ | diffed_comp := NBStrongComponent.SLICED_COMPONENT( | |
| 931 | var_cref = lhsKey, | ||
| 932 | var = Slice.SLICE(lhsVarPtr, {}), | ||
| 933 | eqn = Slice.SLICE(eqPtr, {}), | ||
| 934 | status = NBSolve.Status.EXPLICIT | ||
| 935 | ); | ||
| 936 | else | ||
| 937 | ✗ | diffed_comp := NBStrongComponent.SINGLE_COMPONENT( | |
| 938 | var = lhsVarPtr, | ||
| 939 | eqn = eqPtr, | ||
| 940 | status = NBSolve.Status.EXPLICIT | ||
| 941 | ); | ||
| 942 | end if; | ||
| 943 | then diffed_comp; | ||
| 944 | ✗ | case NBEquation.ARRAY_EQUATION() then | |
| 945 | NBStrongComponent.SINGLE_COMPONENT( | ||
| 946 | var = lhsVarPtr, | ||
| 947 | eqn = eqPtr, | ||
| 948 | status = NBSolve.Status.EXPLICIT | ||
| 949 | ); | ||
| 950 | ✗ | case NBEquation.RECORD_EQUATION() then | |
| 951 | NBStrongComponent.SINGLE_COMPONENT( | ||
| 952 | var = lhsVarPtr, | ||
| 953 | eqn = eqPtr, | ||
| 954 | status = NBSolve.Status.EXPLICIT | ||
| 955 | ); | ||
| 956 | else algorithm | ||
| 957 | ✗ | Error.addMessage(Error.INTERNAL_ERROR, {getInstanceName() + " cannot create adjoint strong component for equation " + NBEquation.Equation.toString(eq)}); | |
| 958 | ✗ | then fail(); | |
| 959 | end match; | ||
| 960 | end makeAdjointComponentFromRhs; | ||
| 961 | |||
| 962 | function addEntryToLPAMap | ||
| 963 | input Pointer<Variable> vptr; | ||
| 964 | input UnorderedMap<ComponentRef, ComponentRef> diff_map; | ||
| 965 | input UnorderedMap<ComponentRef, AdjointTermList> loop_product_adjoint_map; | ||
| 966 | protected | ||
| 967 | Option<ComponentRef> mappedSeed; | ||
| 968 | algorithm | ||
| 969 | 9 | mappedSeed := UnorderedMap.get(BVariable.getVarName(vptr), diff_map); | |
| 970 |
2/4✗ Branch 0 not taken.
✓ Branch 1 taken 9 times.
✓ Branch 2 taken 9 times.
✗ Branch 3 not taken.
|
9 | if isSome(mappedSeed) then |
| 971 | 9 | UnorderedMap.tryAdd(Util.getOption(mappedSeed), {}, loop_product_adjoint_map); | |
| 972 | end if; | ||
| 973 | end addEntryToLPAMap; | ||
| 974 | |||
| 975 | // Resolve base variables that were actually mapped to tmp pDER vars. | ||
| 976 | // This avoids relying on splitOnTrue output ordering semantics. | ||
| 977 | function getBaseTmpVarCandidates | ||
| 978 | input list<NBVariable.VariablePointer> partialVars; | ||
| 979 | input list<NBVariable.VariablePointer> tmpPDerVars; | ||
| 980 | input UnorderedMap<ComponentRef, ComponentRef> diff_map; | ||
| 981 | output list<NBVariable.VariablePointer> baseTmpVars = {}; | ||
| 982 | protected | ||
| 983 | UnorderedSet<ComponentRef> tmpPDerSet; | ||
| 984 | ComponentRef baseCref; | ||
| 985 | Option<ComponentRef> o_mapped; | ||
| 986 | algorithm | ||
| 987 | 4 | tmpPDerSet := UnorderedSet.new(ComponentRef.hash, ComponentRef.isEqual, Util.nextPrime(listLength(tmpPDerVars))); | |
| 988 | |||
| 989 |
2/2✓ Branch 0 taken 11 times.
✓ Branch 1 taken 4 times.
|
15 | for v in tmpPDerVars loop |
| 990 | 11 | UnorderedSet.add(BVariable.getVarName(v), tmpPDerSet); | |
| 991 | end for; | ||
| 992 | |||
| 993 |
2/2✓ Branch 0 taken 22 times.
✓ Branch 1 taken 4 times.
|
26 | for v in partialVars loop |
| 994 | 22 | baseCref := BVariable.getVarName(v); | |
| 995 | 22 | o_mapped := UnorderedMap.get(baseCref, diff_map); | |
| 996 |
5/6✗ Branch 0 not taken.
✓ Branch 1 taken 22 times.
✓ Branch 2 taken 21 times.
✓ Branch 3 taken 1 time.
✓ Branch 6 taken 11 times.
✓ Branch 7 taken 10 times.
|
22 | if isSome(o_mapped) and UnorderedSet.contains(Util.getOption(o_mapped), tmpPDerSet) then |
| 997 | baseTmpVars := v :: baseTmpVars; | ||
| 998 | end if; | ||
| 999 | end for; | ||
| 1000 | |||
| 1001 | 4 | baseTmpVars := listReverse(baseTmpVars); | |
| 1002 | end getBaseTmpVarCandidates; | ||
| 1003 | |||
| 1004 | // Build a filtered diff map for a given variable list. | ||
| 1005 | // For each variable pointer v in 'vars', if there exists a mapping | ||
| 1006 | // base = BVariable.getVarName(v) -> mapped in 'globalDiffMap' | ||
| 1007 | // then add (base -> mapped) to the returned map. | ||
| 1008 | function populateDiffMap | ||
| 1009 | input list<NBVariable.VariablePointer> vars; | ||
| 1010 | input UnorderedMap<ComponentRef, ComponentRef> globalDiffMap; | ||
| 1011 | output UnorderedMap<ComponentRef, ComponentRef> outMap; | ||
| 1012 | protected | ||
| 1013 | ComponentRef baseCref; | ||
| 1014 | Option<ComponentRef> o_mappedCref; | ||
| 1015 | algorithm | ||
| 1016 | 2 | outMap := UnorderedMap.new<ComponentRef>( | |
| 1017 | ComponentRef.hash, ComponentRef.isEqual, Util.nextPrime(listLength(vars)) | ||
| 1018 | ); | ||
| 1019 | |||
| 1020 |
2/2✓ Branch 0 taken 9 times.
✓ Branch 1 taken 2 times.
|
11 | for vp in vars loop |
| 1021 | 9 | baseCref := BVariable.getVarName(vp); | |
| 1022 | 9 | o_mappedCref := UnorderedMap.get(baseCref, globalDiffMap); | |
| 1023 |
2/4✗ Branch 0 not taken.
✓ Branch 1 taken 9 times.
✓ Branch 2 taken 9 times.
✗ Branch 3 not taken.
|
9 | if isSome(o_mappedCref) then |
| 1024 | 9 | UnorderedMap.add(baseCref, Util.getOption(o_mappedCref), outMap); | |
| 1025 | end if; | ||
| 1026 | end for; | ||
| 1027 | end populateDiffMap; | ||
| 1028 | |||
| 1029 | function isSupportedAdjointStrongComponent | ||
| 1030 | input StrongComponent comp; | ||
| 1031 | output Boolean ok; | ||
| 1032 | algorithm | ||
| 1033 | ok := match comp | ||
| 1034 | case StrongComponent.SINGLE_COMPONENT() then true; | ||
| 1035 | case StrongComponent.MULTI_COMPONENT() then true; | ||
| 1036 | case StrongComponent.SLICED_COMPONENT() then true; | ||
| 1037 | case StrongComponent.RESIZABLE_COMPONENT() then true; | ||
| 1038 | case StrongComponent.ALGEBRAIC_LOOP() then true; | ||
| 1039 | 18 | case StrongComponent.ALIAS() then isSupportedAdjointStrongComponent(comp.original); | |
| 1040 | else false; | ||
| 1041 | end match; | ||
| 1042 | end isSupportedAdjointStrongComponent; | ||
| 1043 | |||
| 1044 | type AdjointTermList = list<Expression>; | ||
| 1045 | function generateAdjointComponent | ||
| 1046 | "Generate adjoint strong component(s) for a single primal strong component. | ||
| 1047 | Uses a fresh adjoint_map per component and returns the resulting adjoint | ||
| 1048 | component(s) plus any new temporary variables." | ||
| 1049 | input StrongComponent comp; | ||
| 1050 | input UnorderedMap<ComponentRef, ComponentRef> diff_map; | ||
| 1051 | input UnorderedMap<Path, Function> funcMap; | ||
| 1052 | input Boolean scalarized; | ||
| 1053 | input Boolean staticAsContinuous; | ||
| 1054 | input Pointer<Integer> idx; | ||
| 1055 | input String contextName; | ||
| 1056 | input VariablePointers seedCandidates "for algebraic loop x-inputs"; | ||
| 1057 | input list<Pointer<Variable>> tmpVarCandidates "base tmp variables to also include in diff_map_x for algebraic loops"; | ||
| 1058 | output list<StrongComponent> adjointComps = {}; | ||
| 1059 | output list<Pointer<Variable>> newTmpVars = {}; | ||
| 1060 | protected | ||
| 1061 | StrongComponent c_noalias; | ||
| 1062 | UnorderedMap<ComponentRef, AdjointTermList> fresh_adjoint_map; | ||
| 1063 | Differentiate.DifferentiationArguments diffArgs; | ||
| 1064 | Equation eq; | ||
| 1065 | list<Statement> adjStmts; | ||
| 1066 | Pointer<Equation> eqPtr; | ||
| 1067 | list<Slice<VariablePointer>> adjVarSlices; | ||
| 1068 | // SSA helper: accumulator for pDer vars created for SSA temporaries | ||
| 1069 | Pointer<list<Pointer<Variable>>> ssaPDerVarsPtr = Pointer.create({}); | ||
| 1070 | algorithm | ||
| 1071 | 19 | c_noalias := StrongComponent.removeAlias(comp); | |
| 1072 | |||
| 1073 | () := match c_noalias | ||
| 1074 | local | ||
| 1075 | // ALGEBRAIC_LOOP locals | ||
| 1076 | Tearing tearing; | ||
| 1077 | list<VariablePointer> itVarPtrs; | ||
| 1078 | list<Expression> residuals; | ||
| 1079 | list<Pointer<Variable>> lambdaPtrs; | ||
| 1080 | list<ComponentRef> lambdaCrefs; | ||
| 1081 | Integer iRes; | ||
| 1082 | Pointer<Variable> lhsVarPtr; | ||
| 1083 | ComponentRef newC; | ||
| 1084 | UnorderedMap<ComponentRef, ComponentRef> diff_map_y, diff_map_x, diff_map_union; | ||
| 1085 | UnorderedMap<ComponentRef, AdjointTermList> loop_product_adjoint_map; | ||
| 1086 | list<Pointer<Variable>> seedPtrListX; | ||
| 1087 | list<Pointer<Equation>> linResEqnPtrs; | ||
| 1088 | AdjointTermList terms_j, terms_x; | ||
| 1089 | Expression lhs_j, rhs_j, rhs_x; | ||
| 1090 | Pointer<Equation> resid_j; | ||
| 1091 | Option<ComponentRef> o_ySeedCref, o_pDerX; | ||
| 1092 | ComponentRef ySeedCref, baseX, pDerX; | ||
| 1093 | StrongComponent loopComp; | ||
| 1094 | |||
| 1095 | StrongComponent ssaAlg; | ||
| 1096 | list<tuple<ComponentRef, tuple<ComponentRef, Integer>>> replacements = {}; | ||
| 1097 | list<Pointer<Variable>> newVars = {}; | ||
| 1098 | // SSA seed-init locals (used in MULTI_COMPONENT adjoint) | ||
| 1099 | UnorderedSet<ComponentRef> seenCrefs; | ||
| 1100 | ComponentRef origCref, finalSsaCref, pDerOrigCref, pDerSsaCref; | ||
| 1101 | Type vty; | ||
| 1102 | // x_bar algorithm locals (used in ALGEBRAIC_LOOP adjoint) | ||
| 1103 | list<Statement> xbarStmts; | ||
| 1104 | SizeClassification sc_x; | ||
| 1105 | Operator addOp_x; | ||
| 1106 | Expression accRhs; | ||
| 1107 | // this true when its an initial problem? but we are only in the dynamic case | ||
| 1108 | Boolean init = false; | ||
| 1109 | |||
| 1110 | // ===================== ALGEBRAIC_LOOP ===================== | ||
| 1111 | case StrongComponent.ALGEBRAIC_LOOP(strict = tearing) algorithm | ||
| 1112 | // Collect iteration vars and residual equations and turn into residual expressions | ||
| 1113 | 1 | itVarPtrs := Tearing.getIterationVars(tearing); | |
| 1114 |
4/4✓ Branch 1 taken 1 time.
✓ Branch 2 taken 1 time.
✓ Branch 3 taken 1 time.
✓ Branch 4 taken 1 time.
|
2 | residuals := list(Equation.getResidualExp(Pointer.access(e)) for e in Tearing.getResidualEqns(tearing)); |
| 1115 | |||
| 1116 | // Create scalar lambda_i temporaries | ||
| 1117 | // Is it possible to create it as a vector? | ||
| 1118 | lambdaPtrs := {}; | ||
| 1119 | lambdaCrefs := {}; | ||
| 1120 |
1/2✓ Branch 1 taken 1 time.
✗ Branch 2 not taken.
|
2 | for iIdx in 1:listLength(residuals) loop |
| 1121 | 1 | (lhsVarPtr, newC) := BVariable.makeAuxVar(NBVariable.TEMPORARY_STR, Pointer.access(idx) + 1, Type.REAL(), false); | |
| 1122 | 1 | Pointer.update(idx, Pointer.access(idx) + 1); | |
| 1123 | 1 | (newC, lhsVarPtr) := BVariable.makePDerVar(newC, contextName, isTmp = true); | |
| 1124 | 1 | lambdaPtrs := lhsVarPtr :: lambdaPtrs; | |
| 1125 | 1 | lambdaCrefs := newC :: lambdaCrefs; | |
| 1126 | end for; | ||
| 1127 | 1 | lambdaPtrs := listReverse(lambdaPtrs); | |
| 1128 | 1 | lambdaCrefs := listReverse(lambdaCrefs); | |
| 1129 | newTmpVars := lambdaPtrs; | ||
| 1130 | |||
| 1131 | // Build filtered diff maps | ||
| 1132 | 1 | diff_map_y := populateDiffMap(itVarPtrs, diff_map); | |
| 1133 | 1 | seedPtrListX := listAppend(BVariable.VariablePointers.toList(seedCandidates), tmpVarCandidates); | |
| 1134 |
6/6✓ Branch 2 taken 1 time.
✓ Branch 3 taken 8 times.
✓ Branch 4 taken 9 times.
✓ Branch 5 taken 1 time.
✓ Branch 6 taken 8 times.
✓ Branch 7 taken 1 time.
|
10 | seedPtrListX := list(vp for vp guard(not UnorderedMap.contains(BVariable.getVarName(vp), diff_map_y)) in seedPtrListX); |
| 1135 | 1 | diff_map_x := populateDiffMap(seedPtrListX, diff_map); | |
| 1136 | 1 | diff_map_union := UnorderedMap.merge(diff_map_y, diff_map_x, sourceInfo()); | |
| 1137 | |||
| 1138 | // Pre-populate loop_product_adjoint_map | ||
| 1139 | 1 | loop_product_adjoint_map := UnorderedMap.new<AdjointTermList>(ComponentRef.hash, ComponentRef.isEqual, listLength(itVarPtrs) + listLength(seedPtrListX)); | |
| 1140 |
2/2✓ Branch 1 taken 1 time.
✓ Branch 2 taken 1 time.
|
2 | for vp in itVarPtrs loop addEntryToLPAMap(vp, diff_map_y, loop_product_adjoint_map); end for; |
| 1141 |
2/2✓ Branch 1 taken 8 times.
✓ Branch 2 taken 1 time.
|
9 | for vp in seedPtrListX loop addEntryToLPAMap(vp, diff_map_x, loop_product_adjoint_map); end for; |
| 1142 | |||
| 1143 | // Accumulate reverse-mode adjoints per residual with seed = lambda_i | ||
| 1144 | iRes := 1; | ||
| 1145 |
2/2✓ Branch 0 taken 1 time.
✓ Branch 1 taken 1 time.
|
2 | for residual_i in residuals loop |
| 1146 |
1/2✓ Branch 1 taken 1 time.
✗ Branch 2 not taken.
|
1 | if iRes > listLength(lambdaCrefs) then break; end if; |
| 1147 | 1 | diffArgs := accumulateAdjointForResidual( | |
| 1148 | residual_i, | ||
| 1149 | Expression.fromCref(listGet(lambdaCrefs, iRes)), | ||
| 1150 | diff_map_union, | ||
| 1151 | funcMap, | ||
| 1152 | scalarized, | ||
| 1153 | loop_product_adjoint_map | ||
| 1154 | ); | ||
| 1155 | // Update loop_product_adjoint_map with new adjoint terms collected from this residual | ||
| 1156 | 1 | loop_product_adjoint_map := Util.getOption(diffArgs.adjoint_map); | |
| 1157 | 1 | iRes := iRes + 1; | |
| 1158 | end for; | ||
| 1159 | |||
| 1160 | // Build linear algebraic loop: sum_i(dr_i/dy_j * lambda_i) = y_bar_j | ||
| 1161 | linResEqnPtrs := {}; | ||
| 1162 |
2/2✓ Branch 0 taken 1 time.
✓ Branch 1 taken 1 time.
|
2 | for vp in itVarPtrs loop |
| 1163 | 1 | o_ySeedCref := UnorderedMap.get(BVariable.getVarName(vp), diff_map_y); | |
| 1164 |
2/4✗ Branch 0 not taken.
✓ Branch 1 taken 1 time.
✓ Branch 2 taken 1 time.
✗ Branch 3 not taken.
|
1 | if isSome(o_ySeedCref) then |
| 1165 | 1 | ySeedCref := Util.getOption(o_ySeedCref); | |
| 1166 | 1 | terms_j := UnorderedMap.getOrDefault(ySeedCref, loop_product_adjoint_map, {}); | |
| 1167 | 1 | lhs_j := buildAdjointRhs(ySeedCref, terms_j); | |
| 1168 | 1 | rhs_j := Expression.fromCref(ySeedCref); | |
| 1169 | 1 | resid_j := Equation.makeAssignment(lhs_j, rhs_j, idx, contextName, | |
| 1170 | NBEquation.Iterator.EMPTY(), NBEquation.EquationAttributes.default(NBEquation.EquationKind.CONTINUOUS, false)); | ||
| 1171 | 1 | linResEqnPtrs := Equation.createResidual(resid_j) :: linResEqnPtrs; | |
| 1172 | end if; | ||
| 1173 | end for; | ||
| 1174 | 1 | linResEqnPtrs := listReverse(linResEqnPtrs); | |
| 1175 | |||
| 1176 |
1/2✓ Branch 0 taken 1 time.
✗ Branch 1 not taken.
|
1 | if not listEmpty(linResEqnPtrs) then |
| 1177 | 1 | loopComp := makeLinearAlgebraicLoop(lambdaPtrs, linResEqnPtrs, NONE(), mixed = false, homotopy = false); | |
| 1178 | adjointComps := loopComp :: adjointComps; | ||
| 1179 | end if; | ||
| 1180 | |||
| 1181 | // Build x_bar = -lambda^T * (dr/dx) as a single algorithm component | ||
| 1182 | xbarStmts := {}; | ||
| 1183 |
2/2✓ Branch 0 taken 8 times.
✓ Branch 1 taken 1 time.
|
9 | for seedVarPtrX in seedPtrListX loop |
| 1184 | 8 | baseX := BVariable.getVarName(seedVarPtrX); | |
| 1185 | 8 | o_pDerX := UnorderedMap.get(baseX, diff_map_x); | |
| 1186 |
2/4✗ Branch 0 not taken.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 8 times.
✗ Branch 3 not taken.
|
8 | if isSome(o_pDerX) then |
| 1187 | 8 | pDerX := Util.getOption(o_pDerX); | |
| 1188 | 8 | terms_x := UnorderedMap.getOrDefault(pDerX, loop_product_adjoint_map, {}); | |
| 1189 |
2/2✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
|
8 | if not listEmpty(terms_x) then |
| 1190 | 2 | rhs_x := Expression.negate(buildAdjointRhs(pDerX, terms_x)); | |
| 1191 | 2 | vty := ComponentRef.getComponentType(pDerX); | |
| 1192 |
1/2✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
|
2 | if Expression.containsCref(rhs_x, pDerX) then |
| 1193 | accRhs := rhs_x; | ||
| 1194 | else | ||
| 1195 | 2 | sc_x := sizeClassificationFromType(vty); | |
| 1196 | 2 | addOp_x := Operator.fromClassification((MathClassification.ADDITION, sc_x), vty); | |
| 1197 | 4 | accRhs := SimplifyExp.simplify(Expression.MULTARY({Expression.fromCref(pDerX), rhs_x}, {}, addOp_x)); | |
| 1198 | end if; | ||
| 1199 | 2 | accRhs := Expression.map(accRhs, Expression.repairOperator); | |
| 1200 | 2 | xbarStmts := Statement.ASSIGNMENT( | |
| 1201 | Expression.fromCref(pDerX), accRhs, vty, DAE.emptyElementSource | ||
| 1202 | ) :: xbarStmts; | ||
| 1203 | end if; | ||
| 1204 | end if; | ||
| 1205 | end for; | ||
| 1206 | 1 | xbarStmts := listReverse(xbarStmts); | |
| 1207 |
1/2✓ Branch 0 taken 1 time.
✗ Branch 1 not taken.
|
1 | if not listEmpty(xbarStmts) then |
| 1208 | 1 | eqPtr := Equation.makeAlgorithm(xbarStmts, init); | |
| 1209 | 1 | Equation.createName(eqPtr, idx, contextName); | |
| 1210 | 1 | adjVarSlices := listReverse(collectAdjointVarSlices(xbarStmts, {})); | |
| 1211 | 1 | adjointComps := StrongComponent.MULTI_COMPONENT( | |
| 1212 | vars = adjVarSlices, | ||
| 1213 | eqn = Slice.SLICE(eqPtr, {}), | ||
| 1214 | status = NBSolve.Status.EXPLICIT | ||
| 1215 | ) :: adjointComps; | ||
| 1216 | end if; | ||
| 1217 | then (); | ||
| 1218 | |||
| 1219 | // ===================== SINGLE_COMPONENT (scalar/array/record equation) ===================== | ||
| 1220 | case StrongComponent.SINGLE_COMPONENT() algorithm | ||
| 1221 | 15 | eq := Pointer.access(c_noalias.eqn); | |
| 1222 | |||
| 1223 | // Build fresh adjoint_map | ||
| 1224 | 15 | fresh_adjoint_map := UnorderedMap.new<AdjointTermList>(ComponentRef.hash, ComponentRef.isEqual, 16); | |
| 1225 |
1/2✓ Branch 0 taken 15 times.
✗ Branch 1 not taken.
|
30 | diffArgs := Differentiate.DIFFERENTIATION_ARGUMENTS( |
| 1226 | diffCref = ComponentRef.EMPTY(), | ||
| 1227 | new_vars = {}, | ||
| 1228 | diff_map = SOME(diff_map), | ||
| 1229 | diffType = DifferentiationType.JACOBIAN, | ||
| 1230 | funcMap = funcMap, | ||
| 1231 | scalarized = scalarized, | ||
| 1232 | adjoint_map = SOME(fresh_adjoint_map), | ||
| 1233 | current_grad = Expression.EMPTY(Type.REAL()), | ||
| 1234 | collectAdjoints = true | ||
| 1235 | ); | ||
| 1236 | |||
| 1237 | 15 | (diffArgs, adjStmts) := Differentiate.differentiateEquationAdjoint(eq, diffArgs); | |
| 1238 | |||
| 1239 |
2/2✓ Branch 0 taken 14 times.
✓ Branch 1 taken 1 time.
|
15 | if not listEmpty(adjStmts) then |
| 1240 | 14 | eqPtr := Equation.makeAlgorithm(adjStmts, init); | |
| 1241 | 14 | Equation.createName(eqPtr, idx, contextName); | |
| 1242 | |||
| 1243 | // Collect output variables from adjoint statements (handles FOR and IF nesting) | ||
| 1244 | 14 | adjVarSlices := listReverse(collectAdjointVarSlices(adjStmts, {})); | |
| 1245 | |||
| 1246 | 14 | adjointComps := {StrongComponent.MULTI_COMPONENT( | |
| 1247 | vars = adjVarSlices, | ||
| 1248 | eqn = Slice.SLICE(eqPtr, {}), | ||
| 1249 | status = NBSolve.Status.EXPLICIT | ||
| 1250 | )}; | ||
| 1251 | end if; | ||
| 1252 | then (); | ||
| 1253 | |||
| 1254 | // ===================== MULTI_COMPONENT (algorithm or if-equation) ===================== | ||
| 1255 | case StrongComponent.MULTI_COMPONENT() algorithm | ||
| 1256 | eq := match Pointer.access(Slice.getT(c_noalias.eqn)) | ||
| 1257 | case Equation.ALGORITHM() algorithm | ||
| 1258 | 1 | (ssaAlg, replacements, newVars) := algorithmToSSA(c_noalias); | |
| 1259 |
1/2✗ Branch 1 not taken.
✓ Branch 2 taken 1 time.
|
1 | if Flags.isSet(Flags.DEBUG_ADJOINT) then |
| 1260 | ✗ | print("SSA algorithm for adjoint of component " + StrongComponent.toString(c_noalias) + ":\n" + StrongComponent.toString(ssaAlg) + "\n"); | |
| 1261 | end if; | ||
| 1262 | |||
| 1263 | // ── Register SSA variables in diff_map ── | ||
| 1264 | // For each new SSA variable, create a pDer companion and add the mapping | ||
| 1265 | // ssaCref -> pDerCref to diff_map so the adjoint differentiation can propagate | ||
| 1266 | // gradients through the SSA rename chain (e.g. x_1 -> pDer.x_1). | ||
| 1267 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 1 time.
|
1 | for ssaVarPtr in newVars loop |
| 1268 | ✗ | makeVarTraverse(ssaVarPtr, contextName, ssaPDerVarsPtr, diff_map, | |
| 1269 | function BVariable.makePDerVar(isTmp = true), staticAsContinuous = staticAsContinuous); | ||
| 1270 | end for; | ||
| 1271 | // Collect the newly created pDer vars as temporaries | ||
| 1272 |
1/2✗ Branch 2 not taken.
✓ Branch 3 taken 1 time.
|
1 | for pDerVarPtr in Pointer.access(ssaPDerVarsPtr) loop |
| 1273 | newTmpVars := pDerVarPtr :: newTmpVars; | ||
| 1274 | end for; | ||
| 1275 | |||
| 1276 | then match ssaAlg | ||
| 1277 | 1 | case StrongComponent.MULTI_COMPONENT() then Pointer.access(Slice.getT(ssaAlg.eqn)); | |
| 1278 | ✗ | else Pointer.access(Slice.getT(c_noalias.eqn)); | |
| 1279 | end match; | ||
| 1280 | else algorithm | ||
| 1281 | ✗ | then Pointer.access(Slice.getT(c_noalias.eqn)); | |
| 1282 | end match; | ||
| 1283 | |||
| 1284 | // Build fresh adjoint_map | ||
| 1285 | 1 | fresh_adjoint_map := UnorderedMap.new<AdjointTermList>(ComponentRef.hash, ComponentRef.isEqual, 16); | |
| 1286 |
1/2✓ Branch 0 taken 1 time.
✗ Branch 1 not taken.
|
2 | diffArgs := Differentiate.DIFFERENTIATION_ARGUMENTS( |
| 1287 | diffCref = ComponentRef.EMPTY(), | ||
| 1288 | new_vars = {}, | ||
| 1289 | diff_map = SOME(diff_map), | ||
| 1290 | diffType = DifferentiationType.JACOBIAN, | ||
| 1291 | funcMap = funcMap, | ||
| 1292 | scalarized = scalarized, | ||
| 1293 | adjoint_map = SOME(fresh_adjoint_map), | ||
| 1294 | current_grad = Expression.EMPTY(Type.REAL()), | ||
| 1295 | collectAdjoints = true | ||
| 1296 | ); | ||
| 1297 | |||
| 1298 | 1 | (diffArgs, adjStmts) := Differentiate.differentiateEquationAdjoint(eq, diffArgs); | |
| 1299 | |||
| 1300 | // TODO: Check if it works as intended and make a test case | ||
| 1301 | // ── Prepend seed-initialization statements for the final SSA variable of each | ||
| 1302 | // multi-assigned original variable ── | ||
| 1303 | // The last SSA rename (x_N) represents the final value of x after the algorithm. | ||
| 1304 | // Before the adjoint reverse sweep we must: | ||
| 1305 | // 1. seed pDer.x_N := pDer.x (transfer the incoming gradient for x) | ||
| 1306 | // We iterate replacements in REVERSE line order so the FIRST entry we see for | ||
| 1307 | // each base variable IS its final SSA rename. A local seen-set avoids re-seeding | ||
| 1308 | // non-final renames. | ||
| 1309 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 1 time.
|
1 | if not listEmpty(newVars) then |
| 1310 | ✗ | seenCrefs := UnorderedSet.new(ComponentRef.hash, ComponentRef.isEqual, 4); | |
| 1311 | ✗ | for replacement in listReverse(replacements) loop | |
| 1312 | ✗ | (origCref, (finalSsaCref, _)) := replacement; | |
| 1313 | ✗ | if not UnorderedSet.contains(origCref, seenCrefs) then | |
| 1314 | ✗ | UnorderedSet.add(origCref, seenCrefs); | |
| 1315 | ✗ | if UnorderedMap.contains(origCref, diff_map) and | |
| 1316 | UnorderedMap.contains(finalSsaCref, diff_map) then | ||
| 1317 | ✗ | pDerOrigCref := UnorderedMap.getOrFail(origCref, diff_map); | |
| 1318 | ✗ | pDerSsaCref := UnorderedMap.getOrFail(finalSsaCref, diff_map); | |
| 1319 | ✗ | vty := ComponentRef.getSubscriptedType(pDerSsaCref, true); | |
| 1320 | ✗ | adjStmts := Statement.ASSIGNMENT( | |
| 1321 | Expression.fromCref(pDerSsaCref), | ||
| 1322 | Expression.fromCref(pDerOrigCref), | ||
| 1323 | vty, DAE.emptyElementSource) :: adjStmts; | ||
| 1324 | end if; | ||
| 1325 | end if; | ||
| 1326 | end for; | ||
| 1327 | end if; | ||
| 1328 | |||
| 1329 |
1/2✓ Branch 0 taken 1 time.
✗ Branch 1 not taken.
|
1 | if not listEmpty(adjStmts) then |
| 1330 | 1 | eqPtr := Equation.makeAlgorithm(adjStmts, init); | |
| 1331 | 1 | Equation.createName(eqPtr, idx, contextName); | |
| 1332 | |||
| 1333 | // Collect output variables from adjoint statements (handles FOR and IF nesting) | ||
| 1334 | 1 | adjVarSlices := listReverse(collectAdjointVarSlices(adjStmts, {})); | |
| 1335 | |||
| 1336 | 1 | adjointComps := {StrongComponent.MULTI_COMPONENT( | |
| 1337 | vars = adjVarSlices, | ||
| 1338 | eqn = Slice.SLICE(eqPtr, {}), | ||
| 1339 | status = NBSolve.Status.EXPLICIT | ||
| 1340 | )}; | ||
| 1341 | end if; | ||
| 1342 | then (); | ||
| 1343 | |||
| 1344 | // ===================== ForComponent: SLICED / RESIZABLE / GENERIC ===================== | ||
| 1345 | case StrongComponent.SLICED_COMPONENT() algorithm | ||
| 1346 | 1 | eq := Pointer.access(Slice.getT(c_noalias.eqn)); | |
| 1347 | 1 | adjointComps := generateAdjointForComponent(eq, c_noalias, diff_map, funcMap, scalarized, init, idx, contextName); | |
| 1348 | then (); | ||
| 1349 | |||
| 1350 | case StrongComponent.RESIZABLE_COMPONENT() algorithm | ||
| 1351 | 1 | eq := Pointer.access(Slice.getT(c_noalias.eqn)); | |
| 1352 | 1 | adjointComps := generateAdjointForComponent(eq, c_noalias, diff_map, funcMap, scalarized, init, idx, contextName); | |
| 1353 | then (); | ||
| 1354 | |||
| 1355 | case StrongComponent.GENERIC_COMPONENT() algorithm | ||
| 1356 | ✗ | eq := Pointer.access(Slice.getT(c_noalias.eqn)); | |
| 1357 | ✗ | adjointComps := generateAdjointForComponent(eq, c_noalias, diff_map, funcMap, scalarized, init, idx, contextName); | |
| 1358 | then (); | ||
| 1359 | |||
| 1360 | else algorithm | ||
| 1361 | ✗ | Error.addMessage(Error.INTERNAL_ERROR, {getInstanceName() + " unsupported component type: " + StrongComponent.toString(c_noalias)}); | |
| 1362 | then (); | ||
| 1363 | end match; | ||
| 1364 | end generateAdjointComponent; | ||
| 1365 | |||
| 1366 | function generateAdjointForComponent | ||
| 1367 | "Handle SLICED/RESIZABLE/GENERIC components that wrap for-equations. | ||
| 1368 | Extracts the body equations, differentiates them, wraps in a for-algorithm." | ||
| 1369 | input Equation eq; | ||
| 1370 | input StrongComponent originalComp; | ||
| 1371 | input UnorderedMap<ComponentRef, ComponentRef> diff_map; | ||
| 1372 | input UnorderedMap<Path, Function> funcMap; | ||
| 1373 | input Boolean scalarized; | ||
| 1374 | input Boolean init; | ||
| 1375 | input Pointer<Integer> idx; | ||
| 1376 | input String contextName; | ||
| 1377 | output list<StrongComponent> adjointComps = {}; | ||
| 1378 | protected | ||
| 1379 | UnorderedMap<ComponentRef, AdjointTermList> fresh_adjoint_map; | ||
| 1380 | Differentiate.DifferentiationArguments diffArgs; | ||
| 1381 | list<Statement> adjStmts; | ||
| 1382 | Pointer<Equation> eqPtr; | ||
| 1383 | list<Slice<VariablePointer>> adjVarSlices; | ||
| 1384 | ComponentRef adjVarCref; | ||
| 1385 | algorithm | ||
| 1386 | // Build fresh adjoint_map and diff arguments | ||
| 1387 | 2 | fresh_adjoint_map := UnorderedMap.new<AdjointTermList>(ComponentRef.hash, ComponentRef.isEqual, 16); | |
| 1388 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
4 | diffArgs := Differentiate.DIFFERENTIATION_ARGUMENTS( |
| 1389 | diffCref = ComponentRef.EMPTY(), | ||
| 1390 | new_vars = {}, | ||
| 1391 | diff_map = SOME(diff_map), | ||
| 1392 | diffType = DifferentiationType.JACOBIAN, | ||
| 1393 | funcMap = funcMap, | ||
| 1394 | scalarized = scalarized, | ||
| 1395 | adjoint_map = SOME(fresh_adjoint_map), | ||
| 1396 | current_grad = Expression.EMPTY(Type.REAL()), | ||
| 1397 | collectAdjoints = true | ||
| 1398 | ); | ||
| 1399 | |||
| 1400 | // differentiateEquationAdjoint handles FOR_EQUATION (wraps with reversed iterators) | ||
| 1401 | 2 | (diffArgs, adjStmts) := Differentiate.differentiateEquationAdjoint(eq, diffArgs); | |
| 1402 | |||
| 1403 |
2/2✓ Branch 0 taken 1 time.
✓ Branch 1 taken 1 time.
|
2 | if not listEmpty(adjStmts) then |
| 1404 | 1 | eqPtr := Equation.makeAlgorithm(adjStmts, init); | |
| 1405 | 1 | Equation.createName(eqPtr, idx, contextName); | |
| 1406 | |||
| 1407 | // Collect variable slices from statements (handles ASSIGNMENT, FOR, and IF nesting) | ||
| 1408 | 1 | adjVarSlices := listReverse(collectAdjointVarSlices(adjStmts, {})); | |
| 1409 | |||
| 1410 | // Determine the adjoint var cref for the component wrapper | ||
| 1411 | adjVarCref := match originalComp | ||
| 1412 | case StrongComponent.SLICED_COMPONENT() then originalComp.var_cref; | ||
| 1413 | case StrongComponent.RESIZABLE_COMPONENT() then originalComp.var_cref; | ||
| 1414 | case StrongComponent.GENERIC_COMPONENT() then originalComp.var_cref; | ||
| 1415 | else ComponentRef.EMPTY(); | ||
| 1416 | end match; | ||
| 1417 | |||
| 1418 | 1 | adjointComps := {StrongComponent.MULTI_COMPONENT( | |
| 1419 | vars = adjVarSlices, | ||
| 1420 | eqn = Slice.SLICE(eqPtr, {}), | ||
| 1421 | status = NBSolve.Status.EXPLICIT | ||
| 1422 | )}; | ||
| 1423 | end if; | ||
| 1424 | end generateAdjointForComponent; | ||
| 1425 | |||
| 1426 | function collectAdjointVarSlices | ||
| 1427 | "Recursively collect variable pointer slices from adjoint statements. | ||
| 1428 | Handles ASSIGNMENT at any nesting depth inside FOR and IF bodies." | ||
| 1429 | input list<Statement> stmts; | ||
| 1430 | input output list<Slice<VariablePointer>> varSlices; | ||
| 1431 | protected | ||
| 1432 | Pointer<Variable> vPtr; | ||
| 1433 | ComponentRef baseCref; | ||
| 1434 | algorithm | ||
| 1435 |
2/2✓ Branch 0 taken 33 times.
✓ Branch 1 taken 18 times.
|
51 | for s in stmts loop |
| 1436 | () := match s | ||
| 1437 | case Statement.ASSIGNMENT(lhs = Expression.CREF()) algorithm | ||
| 1438 | 32 | baseCref := ComponentRef.stripSubscriptsAll(Expression.toCref(s.lhs)); | |
| 1439 | try | ||
| 1440 | 32 | vPtr := BVariable.getVarPointer(baseCref, sourceInfo()); | |
| 1441 | 32 | varSlices := Slice.SLICE(vPtr, {}) :: varSlices; | |
| 1442 | else | ||
| 1443 | end try; | ||
| 1444 | then (); | ||
| 1445 | case Statement.FOR() algorithm | ||
| 1446 | 1 | varSlices := collectAdjointVarSlices(s.body, varSlices); | |
| 1447 | then (); | ||
| 1448 | case Statement.IF() algorithm | ||
| 1449 | ✗ | for branch in s.branches loop | |
| 1450 | ✗ | varSlices := collectAdjointVarSlices(Util.tuple22(branch), varSlices); | |
| 1451 | end for; | ||
| 1452 | then (); | ||
| 1453 | else (); | ||
| 1454 | end match; | ||
| 1455 | end for; | ||
| 1456 | end collectAdjointVarSlices; | ||
| 1457 | |||
| 1458 | function jacobianSymbolicAdjoint extends Module.jacobianInterface; | ||
| 1459 | protected | ||
| 1460 | list<StrongComponent> comps, primalComps, diffed_comps = {}; | ||
| 1461 | Pointer<list<Pointer<Variable>>> seed_vars_ptr = Pointer.create({}); | ||
| 1462 | Pointer<list<Pointer<Variable>>> pDer_vars_ptr = Pointer.create({}); | ||
| 1463 | UnorderedMap<ComponentRef,ComponentRef> diff_map = UnorderedMap.new<ComponentRef>(ComponentRef.hash, ComponentRef.isEqual); | ||
| 1464 | Pointer<Integer> idx = Pointer.create(0); | ||
| 1465 | |||
| 1466 | list<Pointer<Variable>> all_vars, unknown_vars, aux_vars, alias_vars, depend_vars, res_vars, tmp_vars, seed_vars, old_res_vars, baseTmpVarCandidates; | ||
| 1467 | BVariable.VarData varDataJac; | ||
| 1468 | |||
| 1469 | VariablePointers adjacencyVars; | ||
| 1470 | Adjacency.Matrix fullLocal, sparsity; | ||
| 1471 | UnorderedSet<ComponentRef> seed_set = UnorderedSet.new(ComponentRef.hash, ComponentRef.isEqual); | ||
| 1472 | UnorderedSet<ComponentRef> pder_set = UnorderedSet.new(ComponentRef.hash, ComponentRef.isEqual); | ||
| 1473 | |||
| 1474 | String newName; | ||
| 1475 | |||
| 1476 | BVariable.checkVar func = getTmpFilterFunction(jacType); | ||
| 1477 | |||
| 1478 | // Per-component adjoint generation | ||
| 1479 | list<StrongComponent> compAdjComps; | ||
| 1480 | list<Pointer<Variable>> compNewVars; | ||
| 1481 | algorithm | ||
| 1482 | 4 | newName := name + "_ADJ"; | |
| 1483 |
2/4✗ Branch 0 not taken.
✓ Branch 1 taken 4 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 4 times.
|
4 | if isSome(strongComponents) then |
| 1484 |
6/6✓ Branch 3 taken 1 time.
✓ Branch 4 taken 19 times.
✓ Branch 5 taken 20 times.
✓ Branch 6 taken 4 times.
✓ Branch 7 taken 19 times.
✓ Branch 8 taken 4 times.
|
47 | comps := list(comp for comp guard(not StrongComponent.isDiscrete(comp)) in Util.getOption(strongComponents)); |
| 1485 | primalComps := comps; | ||
| 1486 | // only allow currently implemented adjoint-capable components | ||
| 1487 |
2/2✓ Branch 0 taken 19 times.
✓ Branch 1 taken 4 times.
|
23 | for c in comps loop |
| 1488 |
1/2✗ Branch 1 not taken.
✓ Branch 2 taken 19 times.
|
19 | if not isSupportedAdjointStrongComponent(c) then |
| 1489 | ✗ | Error.addMessage(Error.INTERNAL_ERROR, { | |
| 1490 | getInstanceName() + " only supports SINGLE_COMPONENT, MULTI_COMPONENT, SLICED_COMPONENT, RESIZABLE_COMPONENT and ALGEBRAIC_LOOP in symbolic adjoint jacobian generation!" | ||
| 1491 | }); | ||
| 1492 | ✗ | fail(); | |
| 1493 | end if; | ||
| 1494 |
1/2✗ Branch 1 not taken.
✓ Branch 2 taken 19 times.
|
19 | if Flags.isSet(Flags.DEBUG_ADJOINT) then |
| 1495 | ✗ | print("Primal component: " + StrongComponent.toString(c) + "\n"); | |
| 1496 | end if; | ||
| 1497 | end for; | ||
| 1498 | else | ||
| 1499 | ✗ | Error.addMessage(Error.INTERNAL_ERROR, {getInstanceName() + " failed because no strong components were given!"}); | |
| 1500 | ✗ | fail(); | |
| 1501 | end if; | ||
| 1502 | |||
| 1503 |
1/2✗ Branch 1 not taken.
✓ Branch 2 taken 4 times.
|
4 | if Flags.isSet(Flags.DEBUG_ADJOINT) then |
| 1504 | ✗ | print("Seed candidates before pDer creation:\n" + BVariable.VariablePointers.toString(seedCandidates, "Seed Candidates") + "\n"); | |
| 1505 | ✗ | print("Partial candidates before pDer creation:\n" + BVariable.VariablePointers.toString(partialCandidates, "Partial Candidates") + "\n"); | |
| 1506 | end if; | ||
| 1507 | |||
| 1508 | // create seed vars | ||
| 1509 |
2/2✓ Branch 1 taken 10 times.
✓ Branch 2 taken 4 times.
|
14 | for v in VariablePointers.toList(seedCandidates) loop |
| 1510 | 10 | makeVarTraverse(v, newName, pDer_vars_ptr, diff_map, function BVariable.makePDerVar(isTmp = false), staticAsContinuous = staticAsContinuous); | |
| 1511 | |||
| 1512 |
1/2✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
|
10 | if BVariable.isContinuous(v, staticAsContinuous) then |
| 1513 | 10 | UnorderedSet.add(BVariable.getVarName(v), seed_set); | |
| 1514 | end if; | ||
| 1515 | end for; | ||
| 1516 | 4 | res_vars := listReverse(Pointer.access(pDer_vars_ptr)); | |
| 1517 | |||
| 1518 | // create pDer vars (also filters out discrete vars) | ||
| 1519 | 4 | (old_res_vars, tmp_vars) := List.splitOnTrue(VariablePointers.toList(partialCandidates), func); | |
| 1520 |
1/2✓ Branch 0 taken 4 times.
✗ Branch 1 not taken.
|
8 | (tmp_vars, _) := List.splitOnTrue(tmp_vars, function BVariable.isContinuous(staticAsContinuous = staticAsContinuous)); |
| 1521 | |||
| 1522 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 4 times.
|
14 | for v in old_res_vars loop |
| 1523 | 10 | UnorderedSet.add(BVariable.getVarName(v), pder_set); | |
| 1524 | end for; | ||
| 1525 | |||
| 1526 |
2/2✓ Branch 1 taken 10 times.
✓ Branch 2 taken 4 times.
|
14 | for v in old_res_vars loop makeVarTraverse(v, newName, seed_vars_ptr, diff_map, BVariable.makeSeedVar, staticAsContinuous = staticAsContinuous); end for; |
| 1527 | 4 | seed_vars := listReverse(Pointer.access(seed_vars_ptr)); | |
| 1528 | |||
| 1529 |
1/2✗ Branch 1 not taken.
✓ Branch 2 taken 4 times.
|
4 | if Flags.isSet(Flags.DEBUG_ADJOINT) then |
| 1530 | ✗ | print("seed vars after seed creation:\n" + BVariable.VariablePointers.toString(VariablePointers.fromList(seed_vars), "Seed Vars") + "\n"); | |
| 1531 | ✗ | print("res vars after pDer creation:\n" + BVariable.VariablePointers.toString(VariablePointers.fromList(res_vars), "Res Vars") + "\n"); | |
| 1532 | ✗ | print("tmp vars after pDer creation:\n" + BVariable.VariablePointers.toString(VariablePointers.fromList(tmp_vars), "Tmp Vars") + "\n"); | |
| 1533 | end if; | ||
| 1534 | |||
| 1535 | 4 | pDer_vars_ptr := Pointer.create({}); | |
| 1536 |
2/2✓ Branch 3 taken 11 times.
✓ Branch 4 taken 4 times.
|
15 | for v in tmp_vars loop makeVarTraverse(v, newName, pDer_vars_ptr, diff_map, function BVariable.makePDerVar(isTmp = true), staticAsContinuous = staticAsContinuous); end for; |
| 1537 | 4 | tmp_vars := Pointer.access(pDer_vars_ptr); | |
| 1538 | 4 | baseTmpVarCandidates := getBaseTmpVarCandidates(VariablePointers.toList(partialCandidates), tmp_vars, diff_map); | |
| 1539 | |||
| 1540 |
1/2✗ Branch 1 not taken.
✓ Branch 2 taken 4 times.
|
4 | if Flags.isSet(Flags.DEBUG_ADJOINT) then |
| 1541 | ✗ | print("Diff map before component generation:\n" + diffMapToString(diff_map) + "\n"); | |
| 1542 | end if; | ||
| 1543 | |||
| 1544 | // ===================== Sequential adjoint component generation ===================== | ||
| 1545 | // Process each primal component in reverse order (LIFO), generate adjoint component(s), | ||
| 1546 | // and prepend to the unified list. | ||
| 1547 |
2/2✓ Branch 0 taken 19 times.
✓ Branch 1 taken 4 times.
|
23 | for comp in primalComps loop |
| 1548 | 19 | (compAdjComps, compNewVars) := generateAdjointComponent( | |
| 1549 | comp, diff_map, funcMap, seedCandidates.scalarized, staticAsContinuous, idx, newName, seedCandidates, baseTmpVarCandidates); | ||
| 1550 | |||
| 1551 | // Prepend adjoint components (already in correct order from generateAdjointComponent) | ||
| 1552 | // only more than one if the original component was an algebraic loop | ||
| 1553 |
2/2✓ Branch 1 taken 18 times.
✓ Branch 2 taken 19 times.
|
37 | for ac in compAdjComps loop |
| 1554 | diffed_comps := ac :: diffed_comps; | ||
| 1555 | end for; | ||
| 1556 | |||
| 1557 | // Collect any new temporary variables (e.g. lambda vars from algebraic loops) | ||
| 1558 |
2/2✓ Branch 0 taken 1 time.
✓ Branch 1 taken 19 times.
|
20 | for v in compNewVars loop |
| 1559 | 1 | tmp_vars := v :: tmp_vars; | |
| 1560 | end for; | ||
| 1561 | |||
| 1562 |
1/2✗ Branch 1 not taken.
✓ Branch 2 taken 19 times.
|
19 | if Flags.isSet(Flags.DEBUG_ADJOINT) then |
| 1563 | ✗ | for ac in compAdjComps loop | |
| 1564 | ✗ | print("[adjoint] generated component: " + StrongComponent.toString(ac) + "\n"); | |
| 1565 | end for; | ||
| 1566 | end if; | ||
| 1567 | end for; | ||
| 1568 | // diffed_comps is now in LIFO order (correct for adjoint execution) | ||
| 1569 | |||
| 1570 |
1/2✗ Branch 1 not taken.
✓ Branch 2 taken 4 times.
|
4 | if Flags.isSet(Flags.DEBUG_ADJOINT) then |
| 1571 | ✗ | print("Final list of differentiated components:\n"); | |
| 1572 | ✗ | for comp in diffed_comps loop | |
| 1573 | ✗ | print(StrongComponent.toString(comp) + "\n"); | |
| 1574 | end for; | ||
| 1575 | end if; | ||
| 1576 | |||
| 1577 | // collect var data (most of this can be removed) | ||
| 1578 | 4 | unknown_vars := listAppend(res_vars, tmp_vars); | |
| 1579 | all_vars := unknown_vars; // add other vars later on | ||
| 1580 | |||
| 1581 | 4 | seed_vars := Pointer.access(seed_vars_ptr); | |
| 1582 | aux_vars := seed_vars; // add other auxiliaries later on. TODO: Need to add the SSA vars and the lambda vars from algebraic loops as auxiliaries? | ||
| 1583 | alias_vars := {}; | ||
| 1584 | depend_vars := {}; | ||
| 1585 | |||
| 1586 | 4 | varDataJac := BVariable.VAR_DATA_JAC( | |
| 1587 | variables = VariablePointers.fromList(all_vars), | ||
| 1588 | unknowns = VariablePointers.fromList(unknown_vars), | ||
| 1589 | auxiliaries = VariablePointers.fromList(aux_vars), | ||
| 1590 | aliasVars = VariablePointers.fromList(alias_vars), | ||
| 1591 | diffVars = partialCandidates, | ||
| 1592 | dependencies = VariablePointers.fromList(depend_vars), | ||
| 1593 | resultVars = VariablePointers.fromList(res_vars), | ||
| 1594 | tmpVars = VariablePointers.fromList(tmp_vars), | ||
| 1595 | seedVars = VariablePointers.fromList(seed_vars) | ||
| 1596 | ); | ||
| 1597 | |||
| 1598 | 4 | adjacencyVars := VariablePointers.clone(seedCandidates); | |
| 1599 | // tmp_vars are diffed so use the undiffed ones for adjacency (but does this adjacency approach even work because the adjoint adds tmp vars to the system which are not part of the original system?) | ||
| 1600 | 4 | adjacencyVars := VariablePointers.addList(baseTmpVarCandidates, adjacencyVars); | |
| 1601 |
1/2✓ Branch 0 taken 4 times.
✗ Branch 1 not taken.
|
4 | if jacType == JacobianType.ODE then |
| 1602 | 4 | adjacencyVars := VariablePointers.addList(VariablePointers.toList(partialCandidates), adjacencyVars); | |
| 1603 | end if; | ||
| 1604 |
4/4✓ Branch 0 taken 19 times.
✓ Branch 1 taken 4 times.
✓ Branch 2 taken 19 times.
✓ Branch 3 taken 4 times.
|
23 | fullLocal := Adjacency.Matrix.createFull(adjacencyVars, |
| 1605 | EquationPointers.fromList(List.flatten(list(StrongComponent.getEquations(comp) for comp in comps)))); | ||
| 1606 | 4 | sparsity := Adjacency.Matrix.fullToSparsity(fullLocal, comps, seed_set, pder_set, diff_map, isAdjoint = true); | |
| 1607 | |||
| 1608 | 4 | jacobian := SOME(Jacobian.JACOBIAN( | |
| 1609 | name = newName, | ||
| 1610 | jacType = jacType, | ||
| 1611 | varData = varDataJac, | ||
| 1612 | comps = listArray(diffed_comps), | ||
| 1613 | sparsity = sparsity, | ||
| 1614 | isAdjoint = true | ||
| 1615 | )); | ||
| 1616 | end jacobianSymbolicAdjoint; | ||
| 1617 | |||
| 1618 | function jacobianNumeric | ||
| 1619 | extends Module.jacobianInterface; | ||
| 1620 | protected | ||
| 1621 | VarData varDataJac; | ||
| 1622 | VariablePointers adjacencyVars; | ||
| 1623 | Adjacency.Matrix sparsity, fullLocal; | ||
| 1624 | list<Pointer<Variable>> res_vars, tmp_vars, seed_vars_d, pDer_vars_d; | ||
| 1625 | BVariable.checkVar func = getTmpFilterFunction(jacType); | ||
| 1626 | Pointer<list<Pointer<Variable>>> seed_vars_ptr = Pointer.create({}); | ||
| 1627 | Pointer<list<Pointer<Variable>>> pDer_vars_ptr = Pointer.create({}); | ||
| 1628 | UnorderedMap<ComponentRef,ComponentRef> diff_map = UnorderedMap.new<ComponentRef>(ComponentRef.hash, ComponentRef.isEqual); | ||
| 1629 | |||
| 1630 | UnorderedSet<ComponentRef> seed_set = UnorderedSet.new(ComponentRef.hash, ComponentRef.isEqual); | ||
| 1631 | UnorderedSet<ComponentRef> pder_set = UnorderedSet.new(ComponentRef.hash, ComponentRef.isEqual); | ||
| 1632 | list<StrongComponent> sparsity_comps; | ||
| 1633 | algorithm | ||
| 1634 | 80 | (res_vars, tmp_vars) := List.splitOnTrue(VariablePointers.toList(partialCandidates), func); | |
| 1635 |
2/2✓ Branch 0 taken 78 times.
✓ Branch 1 taken 2 times.
|
158 | (tmp_vars, _) := List.splitOnTrue(tmp_vars, function BVariable.isContinuous(staticAsContinuous = staticAsContinuous)); |
| 1636 | |||
| 1637 | 80 | VariablePointers.mapPtr(seedCandidates, function makeVarTraverse(name = name, vars_ptr = seed_vars_ptr, map = diff_map, makeVar = BVariable.makeSeedVar, staticAsContinuous = staticAsContinuous)); | |
| 1638 | 80 | seed_vars_d := Pointer.access(seed_vars_ptr); | |
| 1639 |
2/2✓ Branch 1 taken 106 times.
✓ Branch 2 taken 80 times.
|
186 | for v in VariablePointers.toList(seedCandidates) loop |
| 1640 |
1/2✓ Branch 1 taken 106 times.
✗ Branch 2 not taken.
|
106 | if BVariable.isContinuous(v, staticAsContinuous) then |
| 1641 | 106 | UnorderedSet.add(BVariable.getVarName(v), seed_set); | |
| 1642 | // Also add base cref so iterator-subscripted deps from for-loop equations | ||
| 1643 | // can match via base fallback in filterSet. | ||
| 1644 | 106 | UnorderedSet.add(ComponentRef.stripSubscriptsAll(BVariable.getVarName(v)), seed_set); | |
| 1645 | end if; | ||
| 1646 | end for; | ||
| 1647 | |||
| 1648 |
2/2✓ Branch 0 taken 106 times.
✓ Branch 1 taken 80 times.
|
186 | for v in res_vars loop |
| 1649 | 106 | UnorderedSet.add(BVariable.getVarName(v), pder_set); | |
| 1650 | 106 | makeVarTraverse(v, name, pDer_vars_ptr, diff_map, function BVariable.makePDerVar(isTmp = false), staticAsContinuous = staticAsContinuous); | |
| 1651 | end for; | ||
| 1652 | 80 | pDer_vars_d := Pointer.access(pDer_vars_ptr); | |
| 1653 | |||
| 1654 | 80 | varDataJac := BVariable.VAR_DATA_JAC( | |
| 1655 | variables = VariablePointers.fromList({}), | ||
| 1656 | unknowns = partialCandidates, | ||
| 1657 | auxiliaries = VariablePointers.fromList(seed_vars_d), | ||
| 1658 | aliasVars = VariablePointers.fromList({}), | ||
| 1659 | diffVars = partialCandidates, | ||
| 1660 | dependencies = VariablePointers.fromList({}), | ||
| 1661 | resultVars = VariablePointers.fromList(pDer_vars_d), | ||
| 1662 | tmpVars = VariablePointers.fromList(tmp_vars), | ||
| 1663 | seedVars = VariablePointers.fromList(seed_vars_d) | ||
| 1664 | ); | ||
| 1665 | |||
| 1666 |
2/4✗ Branch 0 not taken.
✓ Branch 1 taken 80 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 80 times.
|
80 | if isSome(strongComponents) then |
| 1667 | 80 | adjacencyVars := VariablePointers.clone(seedCandidates); | |
| 1668 | 80 | adjacencyVars := VariablePointers.addList(tmp_vars, adjacencyVars); | |
| 1669 |
2/2✓ Branch 0 taken 76 times.
✓ Branch 1 taken 4 times.
|
80 | if jacType == JacobianType.ODE then |
| 1670 | 76 | adjacencyVars := VariablePointers.addList(res_vars, adjacencyVars); | |
| 1671 | end if; | ||
| 1672 | // with resizable arrays the inner equations of algebraic loops are needed for their dependencies | ||
| 1673 | 80 | sparsity_comps := arrayList(Util.getOption(strongComponents)); | |
| 1674 |
2/2✓ Branch 1 taken 8 times.
✓ Branch 2 taken 72 times.
|
80 | if Flags.getConfigBool(Flags.RESIZABLE_ARRAYS) then |
| 1675 |
4/4✓ Branch 0 taken 66 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 66 times.
✓ Branch 3 taken 8 times.
|
74 | sparsity_comps := List.flatten(list(withInnerComps(comp) for comp in sparsity_comps)); |
| 1676 | end if; | ||
| 1677 |
4/4✓ Branch 0 taken 552 times.
✓ Branch 1 taken 80 times.
✓ Branch 2 taken 552 times.
✓ Branch 3 taken 80 times.
|
632 | fullLocal := Adjacency.Matrix.createFull(adjacencyVars, EquationPointers.fromList( |
| 1678 | List.flatten(list(StrongComponent.getEquations(comp) for comp in sparsity_comps)))); | ||
| 1679 | 80 | sparsity := Adjacency.Matrix.fullToSparsity(fullLocal, sparsity_comps, seed_set, pder_set, diff_map); | |
| 1680 | else | ||
| 1681 | ✗ | Error.addMessage(Error.INTERNAL_ERROR,{getInstanceName() + " failed because strong components are missing."}); | |
| 1682 | ✗ | fail(); | |
| 1683 | end if; | ||
| 1684 | |||
| 1685 | 80 | jacobian := SOME(Jacobian.JACOBIAN( | |
| 1686 | name = name, | ||
| 1687 | jacType = jacType, | ||
| 1688 | varData = varDataJac, | ||
| 1689 | comps = listArray({}), | ||
| 1690 | sparsity = sparsity, | ||
| 1691 | isAdjoint = false | ||
| 1692 | )); | ||
| 1693 | end jacobianNumeric; | ||
| 1694 | |||
| 1695 | function jacobianNone | ||
| 1696 | extends Module.jacobianInterface; | ||
| 1697 | algorithm | ||
| 1698 | jacobian := NONE(); | ||
| 1699 | end jacobianNone; | ||
| 1700 | |||
| 1701 | function getTmpFilterFunction | ||
| 1702 | " - ODE filter by state derivative / algebraic | ||
| 1703 | - LS/NLS/DAE filter by residual / inner" | ||
| 1704 | input JacobianType jacType; | ||
| 1705 | output BVariable.checkVar func; | ||
| 1706 | algorithm | ||
| 1707 | func := match jacType | ||
| 1708 | case JacobianType.ODE then BVariable.isStateDerivative; | ||
| 1709 | case JacobianType.DAE then BVariable.isResidual; | ||
| 1710 | case JacobianType.LS then BVariable.isResidual; | ||
| 1711 | case JacobianType.NLS then BVariable.isResidual; | ||
| 1712 | case JacobianType.OPT_LFG then BVariable.isLfgFunction; | ||
| 1713 | case JacobianType.OPT_MRF then BVariable.isMrfFunction; | ||
| 1714 | case JacobianType.OPT_R0 then BVariable.isInitialConstraint; | ||
| 1715 | else algorithm | ||
| 1716 | ✗ | Error.addMessage(Error.INTERNAL_ERROR,{getInstanceName() + " failed because jacobian type is not known: " + jacobianTypeString(jacType)}); | |
| 1717 | ✗ | then fail(); | |
| 1718 | end match; | ||
| 1719 | end getTmpFilterFunction; | ||
| 1720 | |||
| 1721 | function makeVarTraverse | ||
| 1722 | input Pointer<Variable> var_ptr; | ||
| 1723 | input String name; | ||
| 1724 | input Pointer<list<Pointer<Variable>>> vars_ptr; | ||
| 1725 | input UnorderedMap<ComponentRef,ComponentRef> map; | ||
| 1726 | input Func makeVar; | ||
| 1727 | input Boolean staticAsContinuous; | ||
| 1728 | |||
| 1729 | partial function Func | ||
| 1730 | input output ComponentRef cref; | ||
| 1731 | input String name; | ||
| 1732 | output Pointer<Variable> diff_ptr; | ||
| 1733 | end Func; | ||
| 1734 | protected | ||
| 1735 | Variable var = Pointer.access(var_ptr); | ||
| 1736 | ComponentRef diff, parent_name, diff_parent_name; | ||
| 1737 | Pointer<Variable> diff_ptr, parent, diff_parent; | ||
| 1738 | algorithm | ||
| 1739 | // only create seed or pDer var if it is continuous | ||
| 1740 |
1/2✓ Branch 1 taken 2024 times.
✗ Branch 2 not taken.
|
2024 | if BVariable.isContinuous(var_ptr, staticAsContinuous) then |
| 1741 | // make the new differentiated variable itself | ||
| 1742 |
2/2✓ Branch 0 taken 1170 times.
✓ Branch 1 taken 854 times.
|
2024 | (diff, diff_ptr) := makeVar(var.name, name); |
| 1743 | // add $<new>.x variable pointer to the variables | ||
| 1744 | 2024 | Pointer.update(vars_ptr, diff_ptr :: Pointer.access(vars_ptr)); | |
| 1745 | // add x -> $<new>.x to the map for later lookup | ||
| 1746 | 2024 | UnorderedMap.add(var.name, diff, map); | |
| 1747 | // Base-cref fallback for iterator-subscripted deps (x[$i1]) to find their seed in Part D. | ||
| 1748 | // Literal subscripts (x[1] vs x[2]) are independent unknowns: registering a base-cref | ||
| 1749 | // fallback for them would let an UNRELATED literal element (never itself an unknown of | ||
| 1750 | // this Jacobian, e.g. x[1] when only x[2] is) wrongly resolve via NBDifferentiate's | ||
| 1751 | // exact-match-first lookup falling through to this template. Exclude those; literal | ||
| 1752 | // elements that *are* genuine unknowns already get resolved by the exact-match check. | ||
| 1753 |
3/6✓ Branch 1 taken 122 times.
✓ Branch 2 taken 1902 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 122 times.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
|
2024 | if ComponentRef.hasSubscripts(var.name) and not List.all(ComponentRef.subscriptsAllFlat(var.name), Subscript.isLiteral) |
| 1754 | and not UnorderedMap.contains(ComponentRef.stripSubscriptsAll(var.name), map) then | ||
| 1755 | ✗ | UnorderedMap.add(ComponentRef.stripSubscriptsAll(var.name), diff, map); | |
| 1756 | end if; | ||
| 1757 | |||
| 1758 | // differentiate parent and add to map | ||
| 1759 | () := match BVariable.getParent(var_ptr) | ||
| 1760 | case SOME(parent) algorithm | ||
| 1761 | 92 | parent_name := BVariable.getVarName(parent); | |
| 1762 | diff_parent := match UnorderedMap.get(parent_name, map) | ||
| 1763 | 58 | case SOME(diff_parent_name) then BVariable.getVarPointer(diff_parent_name, sourceInfo()); | |
| 1764 | else algorithm | ||
| 1765 |
2/2✓ Branch 0 taken 19 times.
✓ Branch 1 taken 15 times.
|
34 | (diff_parent_name, _) := makeVar(parent_name, name); |
| 1766 | 34 | UnorderedMap.add(parent_name, diff_parent_name, map); | |
| 1767 | 34 | then BVariable.getVarPointer(diff_parent_name, sourceInfo()); | |
| 1768 | end match; | ||
| 1769 | |||
| 1770 | // add the child to the list of children | ||
| 1771 | 92 | BVariable.addRecordChild(diff_parent, diff_ptr); | |
| 1772 | // set the parent of the child | ||
| 1773 | 92 | diff_ptr := BVariable.setParent(diff_ptr, diff_parent); | |
| 1774 | then (); | ||
| 1775 | |||
| 1776 | else (); | ||
| 1777 | end match; | ||
| 1778 | end if; | ||
| 1779 | end makeVarTraverse; | ||
| 1780 | |||
| 1781 | function diffMapToString | ||
| 1782 | input UnorderedMap<ComponentRef, ComponentRef> map; | ||
| 1783 | output String s; | ||
| 1784 | algorithm | ||
| 1785 | ✗ | s := UnorderedMap.toString(map, ComponentRef.toString, ComponentRef.toString, "\n ", " -> "); | |
| 1786 | ✗ | s := "{\n " + s + "\n}"; | |
| 1787 | end diffMapToString; | ||
| 1788 | |||
| 1789 | function makeLinearAlgebraicLoop | ||
| 1790 | input list<NBVariable.VariablePointer> itVarPtrs; // unknowns y (order = columns of A) | ||
| 1791 | input list<Pointer<NBEquation.Equation>> resEqnPtrs; // residuals r_i(y)=0, same order as rows of A | ||
| 1792 | input Option<NBackendDAE> jac = NONE(); // optional analytic Jacobian for A | ||
| 1793 | input Boolean mixed = false; | ||
| 1794 | input Boolean homotopy = false; | ||
| 1795 | output NBStrongComponent comp; | ||
| 1796 | protected | ||
| 1797 | Integer m1 = listLength(itVarPtrs); | ||
| 1798 | Integer m2 = listLength(resEqnPtrs); | ||
| 1799 | list<NBSlice<NBVariable.VariablePointer>> itVars_s; | ||
| 1800 | list<NBSlice<Pointer<NBEquation.Equation>>> res_s; | ||
| 1801 | NBTearing.Tearing tearingSet; | ||
| 1802 | algorithm | ||
| 1803 | // sanity | ||
| 1804 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 1 time.
|
1 | if m1 <> m2 then |
| 1805 | ✗ | Error.addMessage(Error.INTERNAL_ERROR, {"makeLinearAlgebraicLoop: |vars| != |eqns|"}); | |
| 1806 | ✗ | fail(); | |
| 1807 | end if; | ||
| 1808 | |||
| 1809 | // wrap as full slices (keep order) | ||
| 1810 |
4/4✓ Branch 0 taken 1 time.
✓ Branch 1 taken 1 time.
✓ Branch 2 taken 1 time.
✓ Branch 3 taken 1 time.
|
2 | itVars_s := list(NBSlice.SLICE(vp, {}) for vp in itVarPtrs); |
| 1811 |
4/4✓ Branch 0 taken 1 time.
✓ Branch 1 taken 1 time.
✓ Branch 2 taken 1 time.
✓ Branch 3 taken 1 time.
|
2 | res_s := list(NBSlice.SLICE(ep, {}) for ep in resEqnPtrs); |
| 1812 | |||
| 1813 | // strict tearing: no inner equations for a plain linear system | ||
| 1814 | 1 | tearingSet := NBTearing.TEARING_SET( | |
| 1815 | iteration_vars = itVars_s, | ||
| 1816 | residual_eqns = res_s, | ||
| 1817 | innerEquations = listArray({}), | ||
| 1818 | jac = jac | ||
| 1819 | ); | ||
| 1820 | |||
| 1821 | // mark as linear algebraic loop | ||
| 1822 |
2/4✓ Branch 0 taken 1 time.
✗ Branch 1 not taken.
✓ Branch 2 taken 1 time.
✗ Branch 3 not taken.
|
3 | comp := NBStrongComponent.ALGEBRAIC_LOOP( |
| 1823 | idx = -1, | ||
| 1824 | strict = tearingSet, | ||
| 1825 | casual = NONE(), | ||
| 1826 | linear = true, | ||
| 1827 | mixed = mixed, | ||
| 1828 | homotopy = homotopy, | ||
| 1829 | status = NBSolve.Status.IMPLICIT, | ||
| 1830 | implicitlyCreated = false | ||
| 1831 | ); | ||
| 1832 | end makeLinearAlgebraicLoop; | ||
| 1833 | |||
| 1834 | |||
| 1835 | function makeSSAVar | ||
| 1836 | "Creates a fresh SSA variable named 'baseName_idx' that copies all | ||
| 1837 | attributes from the variable referenced by baseCref. | ||
| 1838 | The new variable and its component reference are linked cyclically | ||
| 1839 | via the InstNode VAR_NODE pointer (same pattern as BVariable.makeAuxVar)." | ||
| 1840 | input ComponentRef baseCref "original base cref (no subscripts)"; | ||
| 1841 | input Integer idx "SSA subscript index (1 for x_1, 2 for x_2, ...)"; | ||
| 1842 | output Pointer<Variable> ssaVarPtr; | ||
| 1843 | output ComponentRef ssaCref; | ||
| 1844 | protected | ||
| 1845 | Pointer<Variable> origVarPtr; | ||
| 1846 | Variable origVar; | ||
| 1847 | InstNode newNode; | ||
| 1848 | Type ty; | ||
| 1849 | algorithm | ||
| 1850 | ✗ | origVarPtr := BVariable.getVarPointer(baseCref, sourceInfo()); | |
| 1851 | ✗ | origVar := Pointer.access(origVarPtr); | |
| 1852 | ✗ | ty := ComponentRef.getSubscriptedType(baseCref, false); | |
| 1853 | |||
| 1854 | // Build a fresh VAR_NODE with the SSA name; the variable pointer is | ||
| 1855 | // initially a dummy and is linked to the real variable by makeVarPtr below. | ||
| 1856 | ✗ | newNode := InstNode.VAR_NODE( | |
| 1857 | ComponentRef.firstName(baseCref) + "_" + intString(idx), | ||
| 1858 | PointerWeak.downgrade(Pointer.createImmutable(NBVariable.DUMMY_VARIABLE))); | ||
| 1859 | ✗ | ssaCref := ComponentRef.fromNode(newNode, ty); | |
| 1860 | |||
| 1861 | // Clear any inherited partner pointers (pDer, seed) so that a fresh pDer | ||
| 1862 | // variable is created for this SSA temporary rather than reusing the | ||
| 1863 | // original variable's existing partner. | ||
| 1864 | ✗ | origVar.backendinfo := BackendInfo.BACKEND_INFO( | |
| 1865 | origVar.backendinfo.varKind, | ||
| 1866 | origVar.backendinfo.attributes, | ||
| 1867 | origVar.backendinfo.annotations, | ||
| 1868 | origVar.backendinfo.var_pre, | ||
| 1869 | NONE() /* var_seed */, | ||
| 1870 | NONE() /* var_pder_res */, | ||
| 1871 | NONE() /* var_pder_tmp */, | ||
| 1872 | origVar.backendinfo.var_start, | ||
| 1873 | origVar.backendinfo.parent | ||
| 1874 | ); | ||
| 1875 | |||
| 1876 | // Establish the cyclic Variable <-> InstNode pointer link | ||
| 1877 | ✗ | (ssaVarPtr, ssaCref) := BVariable.makeVarPtr(origVar, ssaCref); | |
| 1878 | end makeSSAVar; | ||
| 1879 | |||
| 1880 | function algorithmToSSA | ||
| 1881 | "Transforms a MULTI_COMPONENT algorithm strong component into SSA | ||
| 1882 | (Static Single Assignment) form. | ||
| 1883 | |||
| 1884 | Variables assigned more than once receive fresh indexed names, | ||
| 1885 | e.g. x -> x_1, x_2, ... RHS reads are updated to use the latest | ||
| 1886 | SSA name of each written variable. | ||
| 1887 | Only ASSIGNMENT statements are expected in the algorithm body. | ||
| 1888 | |||
| 1889 | Each entry (orig_cref, (ssa_cref, line_index)) in `replacements` | ||
| 1890 | records that orig_cref was renamed to ssa_cref at the statement | ||
| 1891 | with 1-based index line_index within the original algorithm." | ||
| 1892 | input StrongComponent comp; | ||
| 1893 | output StrongComponent ssaComp; | ||
| 1894 | output list<tuple<ComponentRef, tuple<ComponentRef, Integer>>> replacements | ||
| 1895 | "original_var -> (ssa_var, line_of_replacement)"; | ||
| 1896 | output list<Pointer<Variable>> newVars | ||
| 1897 | "newly created SSA variable pointers; caller must register them in the variable system"; | ||
| 1898 | protected | ||
| 1899 | Equation eqn; | ||
| 1900 | Algorithm alg; | ||
| 1901 | Statement stmt; | ||
| 1902 | ComponentRef lhsCref, baseCref, ssaCref; | ||
| 1903 | Integer cnt, idx, lineIdx; | ||
| 1904 | Pointer<Variable> ssaVarPtr; | ||
| 1905 | Expression lhsExp, rhsExp; | ||
| 1906 | // Phase 1: how many times is each base cref assigned? | ||
| 1907 | UnorderedMap<ComponentRef, Integer> assignCount = | ||
| 1908 | UnorderedMap.new<Integer>(ComponentRef.hash, ComponentRef.isEqual); | ||
| 1909 | // Phase 2: current per-variable SSA counter | ||
| 1910 | UnorderedMap<ComponentRef, Integer> ssaIdx = | ||
| 1911 | UnorderedMap.new<Integer>(ComponentRef.hash, ComponentRef.isEqual); | ||
| 1912 | // Phase 2: current active SSA expression for each multi-assigned cref | ||
| 1913 | UnorderedMap<ComponentRef, Expression> activeRepl = | ||
| 1914 | UnorderedMap.new<Expression>(ComponentRef.hash, ComponentRef.isEqual); | ||
| 1915 | list<Statement> ssaStmts = {}; | ||
| 1916 | list<tuple<ComponentRef, tuple<ComponentRef, Integer>>> replAcc = {}; | ||
| 1917 | list<Pointer<Variable>> newVarsAcc = {}; | ||
| 1918 | Pointer<Equation> ssaEqnPtr; | ||
| 1919 | algorithm | ||
| 1920 | (ssaComp, replacements, newVars) := match comp | ||
| 1921 | |||
| 1922 | case StrongComponent.MULTI_COMPONENT() algorithm | ||
| 1923 | 1 | eqn := Pointer.access(Slice.getT(comp.eqn)); | |
| 1924 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 1 time.
|
1 | Equation.ALGORITHM(alg = alg) := eqn; |
| 1925 | |||
| 1926 | // ── Phase 1: count how many times each base cref appears on the LHS ── | ||
| 1927 |
2/2✓ Branch 0 taken 3 times.
✓ Branch 1 taken 1 time.
|
4 | for origStmt in alg.statements loop |
| 1928 | () := match origStmt | ||
| 1929 | case Statement.ASSIGNMENT() algorithm | ||
| 1930 | lhsCref := match origStmt.lhs | ||
| 1931 | case Expression.CREF(cref = lhsCref) then lhsCref; | ||
| 1932 | else ComponentRef.EMPTY(); | ||
| 1933 | end match; | ||
| 1934 |
1/2✓ Branch 1 taken 3 times.
✗ Branch 2 not taken.
|
3 | if not ComponentRef.isEmpty(lhsCref) then |
| 1935 | 3 | baseCref := ComponentRef.stripSubscriptsAll(lhsCref); | |
| 1936 | 3 | cnt := UnorderedMap.getOrDefault(baseCref, assignCount, 0); | |
| 1937 | 3 | UnorderedMap.add(baseCref, cnt + 1, assignCount); | |
| 1938 | end if; | ||
| 1939 | then (); | ||
| 1940 | else (); | ||
| 1941 | end match; | ||
| 1942 | end for; | ||
| 1943 | |||
| 1944 | // ── Phase 2: rename multi-assigned variables; substitute RHS reads ── | ||
| 1945 | lineIdx := 1; | ||
| 1946 |
2/2✓ Branch 0 taken 3 times.
✓ Branch 1 taken 1 time.
|
4 | for origStmt in alg.statements loop |
| 1947 | stmt := match origStmt | ||
| 1948 | case Statement.ASSIGNMENT() algorithm | ||
| 1949 | // Substitute every RHS read with its current SSA name | ||
| 1950 | 3 | rhsExp := Expression.map(origStmt.rhs, | |
| 1951 | function Replacements.applySimpleExp(replacements = activeRepl)); | ||
| 1952 | |||
| 1953 | // Check whether the LHS variable needs SSA renaming | ||
| 1954 | 3 | lhsExp := origStmt.lhs; | |
| 1955 | lhsCref := match origStmt.lhs | ||
| 1956 | case Expression.CREF(cref = lhsCref) then lhsCref; | ||
| 1957 | else ComponentRef.EMPTY(); | ||
| 1958 | end match; | ||
| 1959 | |||
| 1960 |
1/2✓ Branch 1 taken 3 times.
✗ Branch 2 not taken.
|
3 | if not ComponentRef.isEmpty(lhsCref) then |
| 1961 | 3 | baseCref := ComponentRef.stripSubscriptsAll(lhsCref); | |
| 1962 |
1/2✗ Branch 1 not taken.
✓ Branch 2 taken 3 times.
|
3 | if UnorderedMap.getOrDefault(baseCref, assignCount, 1) > 1 then |
| 1963 | // Increment the SSA index and create a fresh variable | ||
| 1964 | ✗ | idx := UnorderedMap.getOrDefault(baseCref, ssaIdx, 0) + 1; | |
| 1965 | ✗ | UnorderedMap.add(baseCref, idx, ssaIdx); | |
| 1966 | ✗ | (ssaVarPtr, ssaCref) := makeSSAVar(baseCref, idx); | |
| 1967 | newVarsAcc := ssaVarPtr :: newVarsAcc; | ||
| 1968 | |||
| 1969 | // Re-attach original subscripts to the new SSA cref | ||
| 1970 | ✗ | ssaCref := ComponentRef.copySubscripts(lhsCref, ssaCref); | |
| 1971 | |||
| 1972 | // Update active replacement map (keyed by unsubscripted base cref) | ||
| 1973 | ✗ | UnorderedMap.add(baseCref, | |
| 1974 | Expression.fromCref(ComponentRef.stripSubscriptsAll(ssaCref)), | ||
| 1975 | activeRepl); | ||
| 1976 | |||
| 1977 | // Record: original base cref -> (ssa base cref, 1-based line index) | ||
| 1978 | ✗ | replAcc := (baseCref, | |
| 1979 | (ComponentRef.stripSubscriptsAll(ssaCref), lineIdx)) :: replAcc; | ||
| 1980 | |||
| 1981 | // Replace the LHS with the SSA cref expression | ||
| 1982 | ✗ | lhsExp := Expression.fromCref(ssaCref); | |
| 1983 | end if; | ||
| 1984 | end if; | ||
| 1985 | 3 | then Statement.ASSIGNMENT(lhsExp, rhsExp, origStmt.ty, origStmt.source); | |
| 1986 | |||
| 1987 | else origStmt; | ||
| 1988 | end match; | ||
| 1989 | |||
| 1990 | ssaStmts := stmt :: ssaStmts; | ||
| 1991 | 3 | lineIdx := lineIdx + 1; | |
| 1992 | end for; | ||
| 1993 | |||
| 1994 | // Build a fresh equation pointer with the SSA statement list so the | ||
| 1995 | // original primal equation is left untouched. Is that intended? SSA variables are appended to the component's var list so | ||
| 1996 | // that code generation can declare them as local temporaries. | ||
| 1997 | 1 | alg.statements := listReverse(ssaStmts); | |
| 1998 | eqn := match eqn | ||
| 1999 | 1 | case Equation.ALGORITHM() algorithm eqn.alg := alg; then eqn; | |
| 2000 | else eqn; | ||
| 2001 | end match; | ||
| 2002 | 1 | ssaEqnPtr := Pointer.create(eqn); | |
| 2003 |
2/4✗ Branch 1 not taken.
✓ Branch 2 taken 1 time.
✗ Branch 3 not taken.
✓ Branch 4 taken 1 time.
|
1 | then (StrongComponent.MULTI_COMPONENT( |
| 2004 | vars = listAppend(comp.vars, list(Slice.SLICE(v, {}) for v in listReverse(newVarsAcc))), | ||
| 2005 | eqn = Slice.SLICE(ssaEqnPtr, {}), | ||
| 2006 | status = comp.status | ||
| 2007 | ), listReverse(replAcc), listReverse(newVarsAcc)); | ||
| 2008 | |||
| 2009 | else algorithm | ||
| 2010 | ✗ | Error.addMessage(Error.INTERNAL_ERROR, | |
| 2011 | {getInstanceName() + " expects a MULTI_COMPONENT with an ALGORITHM equation."}); | ||
| 2012 | ✗ | then fail(); | |
| 2013 | |||
| 2014 | end match; | ||
| 2015 | end algorithmToSSA; | ||
| 2016 | |||
| 2017 | annotation(__OpenModelica_Interface="nbackend"); | ||
| 2018 | end NBJacobian; | ||
| 2019 |