casadi::IndexReduction Namespace Reference

Typedefs
typedef struct VariableStruct	Variable
typedef struct EquationStruct	Equation
typedef struct GraphStruct	Graph

Functions
void	graph_add_der (Graph &G, Variable *v)
void	add_variable (Equation *e, Variable *v)
void	graph_add_der (Graph &G, Equation *e, bool add_old=false)
bool	dfs_match_pantelides (Equation *i)
bool	dae_struct_detect_pantelides (Graph &G, std::vector< casadi_int > &var_map, std::vector< casadi_int > &eq_map, casadi_int max_iter)
	Perform Pantelides algorithm for DAE index reduction. More...
void	dae_struct_detect (const std::string &algorithm, const Sparsity &graph, std::vector< casadi_int > &var_map, std::vector< casadi_int > &eq_map, casadi_int max_iter)

Typedef Documentation

Equation

typedef struct EquationStruct casadi::IndexReduction::Equation

Definition at line 1226 of file integration_tools.cpp.

Graph

typedef struct GraphStruct casadi::IndexReduction::Graph

Definition at line 1226 of file integration_tools.cpp.

Variable

typedef struct VariableStruct casadi::IndexReduction::Variable

Definition at line 1226 of file integration_tools.cpp.

Function Documentation

add_variable()

void casadi::IndexReduction::add_variable	(	Equation *	e,
		Variable *	v
	)

Definition at line 478 of file integration_tools.cpp.

                                              {
    auto it = std::find(e->vars0.begin(), e->vars0.end(), v);
    if (it==e->vars0.end()) {
      e->vars0.push_back(v);
      if (!v->deleted) {
        e->vars.push_back(v);
        v->eqs.push_back(e);
      }
    }
  }

Referenced by graph_add_der().

dae_struct_detect()

void casadi::IndexReduction::dae_struct_detect	(	const std::string &	algorithm,
		const Sparsity &	graph,
		std::vector< casadi_int > &	var_map,
		std::vector< casadi_int > &	eq_map,
		casadi_int	max_iter
	)

Definition at line 655 of file integration_tools.cpp.

                         {
    // Input sanitization
    casadi_assert(var_map.size()==graph.size2(),
      "var_map size must match graph columns.");
 
    // Graph structure
    Graph G;
 
    // Allocate space for node representation
    G.E.resize(graph.size1(), nullptr);
    G.V.resize(graph.size2(), nullptr);
    G.nrow_orig = graph.size1();
    G.ncol_orig = graph.size2();
 
    for (auto*& e : G.E) e = new Equation();
    for (auto*& v : G.V) v = new Variable();
 
    // Set positions
    casadi_int i;
    i=0; for (auto* e : G.E) e->i = i++;
    i=0; for (auto* v : G.V) v->i = i++;
 
    // Create edges using incidence sparsity
    const casadi_int* colind = graph.colind();
    const casadi_int* row = graph.row();
 
    // Loop over incidence columns
    for (casadi_int c=0;c<graph.size2();++c) {
      // Loop over nonzeros
      for (casadi_int j=colind[c]; j<colind[c+1]; ++j) {
        casadi_int r = row[j];
        G.V[c]->eqs.push_back(G.E[r]); // Non-monotone
        G.E[r]->vars.push_back(G.V[c]);
      }
    }
 
    // Process var_map
    for (casadi_int i=0;i<var_map.size();++i) {
      Variable* v = G.V[i];
      casadi_int i_der = var_map[i];
      if (i_der>=0) {
        Variable* e = G.V[i_der];
        v->der = e;
        e->der_inv = v;
      }
    }
 
    for (auto* v : G.V) v->eqs0 = v->eqs;
    for (auto* e : G.E) e->vars0 = e->vars;
 
    bool detect = false;
    if (algorithm=="pantelides") {
      detect = dae_struct_detect_pantelides(G, var_map, eq_map, max_iter);
    } else {
      // Clean up
      for (auto* v : G.V) delete v;
      for (auto* e : G.E) delete e;
      casadi_error("Algorithm '" + algorithm + "' not recognized.");
    }
    if (!detect) {
      for (auto* v : G.V) delete v;
      for (auto* e : G.E) delete e;
      casadi_error("Structural detection failed.");
    }
 
    // Prepare outputs
    var_map.resize(G.V.size());
    eq_map.resize(G.E.size());
 
    int k = 0;
    // Populate var_map
    for (auto* v : G.V) {
      if (v->der) {
        var_map[v->i] = v->der->i;
      } else  {
        var_map[v->i] = -1;
      }
      casadi_assert_dev(v->i==k);
      k++;
    }
 
    k = 0;
    // Populate eq_map
    for (auto* e : G.E) {
      if (e->der) {
        eq_map[e->i] = e->der->i;
      } else  {
        eq_map[e->i] = -1;
      }
      casadi_assert_dev(e->i==k);
      k++;
    }
 
    // Clean up
    for (auto* v : G.V) delete v;
    for (auto* e : G.E) delete e;
 
  }

References casadi::Sparsity::colind(), dae_struct_detect_pantelides(), casadi::Sparsity::row(), casadi::Sparsity::size1(), and casadi::Sparsity::size2().

Referenced by casadi::reduce_index_gen().

dae_struct_detect_pantelides()

bool casadi::IndexReduction::dae_struct_detect_pantelides	(	Graph &	G,
		std::vector< casadi_int > &	var_map,
		std::vector< casadi_int > &	eq_map,
		casadi_int	max_iter
	)

The algorithm works purely on structure: no symbolics equations are used.

Parameters

graph	Structural relation between equations (columns) and variables (rows)
var_map	Indicate for each variable where its derivative can be found in the variable list if var_map[i]>=0: var[var_map[i]] == dot(var[i]) Size will increase to accommodate new variables introduced by index reduction
eq_map	Indicate for each equation if it should be differentiated if eq_map[i]>=0: eq[eq_map[i]] == dot(eq[i]) Size will increase over the original number of equations to accommodate extra equations introduced by index reduction

Definition at line 575 of file integration_tools.cpp.

                         {
 
    // Pantelides alg 4.1: Set N'=N do (2)
    casadi_int Np = G.E.size();
    // Pantelides alg 4.1: For k=1 to N' do (3)
    for (casadi_int k=0;k<Np;++k) {
      // Pantelides alg 4.1: Set i=k (3a)
      Equation* i = G.E[k];
 
      casadi_int iter = 0;
      // Pantelides alg 4.1: Repeat (3b)
      while (true) {
        if (iter++ > max_iter) return false;
        // Pantelides alg 4.1: delete all V-nodes with A!=0 and their incident edges
        //       * from the graph (3b-1)
        for (auto* v : G.V) {
          if (v->der) {
            v->deleted = true;
            v->eqs.clear();
            for (auto* e : G.E) {
              auto it = std::find(e->vars.begin(), e->vars.end(), v);
              if (it!=e->vars.end()) e->vars.erase(it);
            }
          }
        }
 
        // Pantelides alg 4.1: Designate all nodes as "uncoloured" (3b-2)
        for (auto* v : G.V) v->visited = false;
        for (auto* e : G.E) e->visited = false;
 
        // Pantelides alg 4.1: Set PATHFOUND= FALSE        (3b-3)
        // Pantelides alg 4.1: AUGMENTPATH (i, PATHFOUND)  (3b-4)
        // Pantelides alg 4.1: Until PATHFOUND  (3c)
        if (dfs_match_pantelides(i)) break; // out of while loop
 
        // Pantelides alg 4.1: If PATHFOUND FALSE then (3b-5)
        casadi_int n = G.V.size();
        // Pantelides alg 4.1: For every coloured V-node j do (3b-5-i)
        for (casadi_int jj=0;jj<n;++jj) {
          Variable* j = G.V[jj];
          if (j->visited && !j->deleted) {
            graph_add_der(G, j);
          }
        }
 
        // Pantelides alg 4.1: For every coloured E-node l do (3b-5-ii)
        n = G.E.size();
        for (casadi_int ll=0;ll<n;++ll) {
          Equation* l = G.E[ll];
          if (l->visited) graph_add_der(G, l, true);
        }
 
        for (casadi_int ll=n;ll<G.E.size();++ll) {
          Equation* l = G.E[ll];
          bool valid = false;
          for (auto* v : l->vars0) {
            if (!v->eqs.empty()) valid = true;
          }
          if (!valid) return false;
        }
 
        // Pantelides alg 4.1: For every coloured V-node j set
        // ASSIGN (A(j))= B(ASSIGN (j))  (3b-5-iii)
        for (auto* j : G.V) {
          if (j->visited && !j->deleted) {
            j->der->assign = j->assign->der;
            j->assign->der->assign = j->der;
          }
        }
 
        // Pantelides alg 4.1: Set i=B(i) (3b-5-iv)
        i = i->der;
 
      }
    }
    return true;
  }

References dfs_match_pantelides(), and graph_add_der().

Referenced by dae_struct_detect().

dfs_match_pantelides()

bool casadi::IndexReduction::dfs_match_pantelides

(

Equation *

i

)

Definition at line 515 of file integration_tools.cpp.

                                         {
 
    // Pantelides alg 3.2: Colour i (1)
    i->visited = true;
 
    // Look for unassigned candidate variables
    // Pantelides alg 3.2: (2) If a V-node j exists such that edge (i-j) exists...
    for (auto* j : i->vars) {
      // Pantelides alg 3.2: ... and ASSIGN (j)= 0 then:
      if (j->assign==nullptr) {
        // Found such variable
        // Assign
        // Pantelides alg 3.2: (2b) Set ASSIGN(j)=i
        j->assign = i; i->assign = j;
 
        // Pantelides alg 3.2: (2a) set PATHFOUND = TRUE, (2c) return
        return true; // End successfully
      }
    }
 
    // Look for assigned candidate variables
    // Pantelides alg 3.2: (3) For every j such that edge (i-j) exists
    for (auto j : i->vars) {
      // Pantelides alg 3.2: ... and j is uncoloured do:
      if (!j->visited) {
        // Pantelides alg 3.2: (3a) Colour j
        j->visited = true;
 
        // Pantelides alg 3.2: (3b) Set k = ASSIGN (j)
        Equation* k = j->assign;
 
        // Pantelides alg 3.2: (3c) AUGMENTPATH (k, PATHFOUND)
        if (dfs_match_pantelides(k)) {
          // Re-assignment
          // Pantelides alg 3.2: (3d-1) Set ASSIGN (j) = i
          j->assign = i; i->assign = j;
          // Pantelides alg 3.2: (3d-2) Return
          return true;
        }
 
      }
    }
 
    // Exhausted all options; return unsuccessfully
    return false;
  }

Referenced by dae_struct_detect_pantelides().

graph_add_der() [1/2]

void casadi::IndexReduction::graph_add_der	(	Graph &	G,
		Equation *	e,
		bool	add_old = false
	)

Definition at line 489 of file integration_tools.cpp.

                                                                {
    // Push new equation to the graph
    G.E.push_back(new Equation());
 
    Equation* e_new = G.E.back();
 
    // Get a position
    e_new->i = G.E.size()-1;
 
    // Add derivative relationship between new and old equation
    e_new->der_inv = e;
    e->der = e_new;
 
    // Loop over all variables in old equation
    for (auto* v : e->vars0) {
      // We depend on the old variable
      if (add_old) {
        add_variable(e_new, v);
      }
      // And we depend on its derivative (create if not present)
      if (!v->der) graph_add_der(G, v);
      add_variable(e_new, v->der);
    }
 
  }

References add_variable(), and graph_add_der().

graph_add_der() [2/2]

void casadi::IndexReduction::graph_add_der	(	Graph &	G,
		Variable *	v
	)

Definition at line 465 of file integration_tools.cpp.

                                            {
    // Push new equation to the graph
    G.V.push_back(new Variable());
    Variable* v_new = G.V.back();
 
    // Get a position
    v_new->i = G.V.size()-1;
 
    // Add derivative relationship between new and old variable
    v_new->der_inv = v;
    v->der = v_new;
  }

Referenced by dae_struct_detect_pantelides(), and graph_add_der().