WebSVN – wimsdev – Diff – /trunk/wims/src/Misc/chemeq/src/chemeq.cc

+// -*- coding: utf-8 -*-
 #include "chemeq.h"
 #include <math.h>
 #include <sstream>
+#include <stdlib.h>
 atome lesatomes[] ={
 {-1, "e"},
 {1, "H"},
 {2, "He"},
+}
 double AtomeListe::weight(fraction mult)const{
   const AtomeListe * al;
   double w=0.0;
   if(Zed!=0 && Zed!=-1 ){ /* cas où ce n'est pas un groupe ou un électron */
     w=mendelweight(symb)*nb*mult.i/mult.d;
-    //debug(); std::cout << std::endl;
+  }
   else if (Zed==0){ /* cas d'un groupe */
     if (group) w+=group->weight(mult*nb);
+  }
   if (suiv) w+=suiv->weight(mult);
   return w;
+}
 void AtomeListe::compte(Compteur &c, fraction mult)const{
   const AtomeListe * al;
   if(Zed!=0 && Zed!=-1 ){ /* cas où ce n'est pas un groupe ou un électron */
     std::string key(symb);
     c[key] +=1.0*nb*mult.i/mult.d;
+  }
   else if (Zed==0){ /* cas d'un groupe */
     if (group) group->compte(c,mult*nb);
+  }
   if (suiv) suiv->compte(c,mult);
+}
 void AtomeListe::numerote(int n){
   if(Zed!=0){ /* cas où ce n'est pas un groupe */
     no = n;
+  }
   else if (Zed==0){ /* cas d'un groupe */
     no = n;
     if (group) group->numerote();
   if (suiv) suiv->numerote(n+1);
+}
 AtomeListe * AtomeListe::triage(AtomeListe * al){
   AtomeListe * al1;
   if(al->Z()!=0){ /* cas où ce n'est pas un groupe */
     if (al->suiv){
       al->suiv = triage(al->suiv);
+    }
     while (al->suiv && al->suiv->Zed!=0 &&
            strcmp(al->symbole(), al->suiv->symbole()) > 1){
 bool operator == (const AtomeListe & a1, const AtomeListe & a2){
   std::stringstream s1, s2;
   a1.printnorm(s1);
   a2.printnorm(s2);
-  return s1 == s2;
+  return s1.str() == s2.str();
+}
 void AtomeListe::printcount(std::ostream & o, const fraction& n, int multiple=1) const{
   if(Zed!=0){ /* cas où ce n'est pas un groupe */
     o << symb;
     o << ':' << n << '*' << multiple*nb;
+  }
   else{ /* cas d'un groupe */
     if (group) group->printcount(o,n,nb);
   if (suiv) {o << ' '; suiv->printcount(o,n,multiple);}
+}
 void AtomeListe::printnorm(std::ostream & o) const{
   if (sqbr) o << "[";
   if(Zed!=0){ /* cas où ce n'est pas un groupe */
     o << symb;
     if (nb!=1) o << nb;
+  }
   else{ /* cas d'un groupe */
     o << "(";
   const AtomeListe * al;
   if (l.sq()) o << "[";
   if(l.Z()>0 || l.Z()<-1){
     o << l.symbole();
+  }
   else if (l.Z()==-1){ // cas de l'électron
     o << "e";
+  }
   else{                // cas des groupes parenthésés
     o << "(";
     if((al=l.groupe())) o << *al;
     o << ")";
+  }
   if((n=l.getmolecularite())>1) o << "_{" << n << "}";
   if((al=l.suivant())) o << *al;
   if (l.sq()) o << "]";
   return o;
+}
-char* moltypeStr[] = { "aq", "g", "s" };
+const char* moltypeStr[] = { "aq", "g", "s" };
 bool operator == (const Molec & m1, const Molec & m2){
   return *(m1.al) == *(m2.al) && m1.ch == m2.ch;
+}
   nb.i *= f.i;
   nb.d *= f.d;
   nb.simplifie();
+}
-bool Molec::printNernst(std::ostream & o, char * prefix){
+bool Molec::printNernst(std::ostream & o, const char * prefix){
   switch(t){
   case sol : return 0;
   case aqueous :
     if (iswater() || iselectron()) return 0;
     o << prefix << "[" << *al;
+}
 fraction Molec::nbelectron()const{
   if (iselectron()) return nb;
   else return fraction(0);
+}
+void Molec::add(fraction f){
+  nb = nb+f;
+}
+void Molec::sub(fraction f){
+  nb = nb-f;
+}
 void Molec::printNombre(std::ostream & o)const{
   if (nb.d==1){
     o << nb.i << "\\,";
+  }
   else {
     o << "\\frac{" << nb.i << "}{" << nb.d << "}\\,";
   }
 }
 std::ostream & operator << (std::ostream & o, const Molec & m){
   if (m.nombre() != 1) m.printNombre(o);
   o << m.liste();
   if (m.charge()){
     o << "^{";
     if(fabs(1.0*m.charge())!=1) o << fabs(1.0*m.charge());
     if(m.charge()>0) o << "+}"; else o << "-}";
   }
   if (m.typage() != aqueous) o << "_{" << moltypeStr[m.typage()] << "}";
   return o;
+}
+int Membre::findMol(const Molec * m){
+  // returns the index of a molecule with the same atomlist if any
+  // else returns -1.
+  int result=-1;
+  for(int i=0; i<size(); i++){
+    if ((*this)[i]->eqMol(m)) result=i;
+  }
+  return result;
+}
+void Membre::addMol(const Molec * m){
+  int i = findMol(m);
+  if (i < 0){
+    push_back(new Molec(*m));
+  } else {
+    (*this)[i]->add(m->nombre());
+  }
+}
+void Membre::addMembre(const Membre * m){
+  for(int i=0; i<m->size(); i++){
+    addMol((*m)[i]);
+  }
+}
+void Membre::eraseNull(){
+  Membre m(*this);
+  clear();
+  for(int i=0; i < m.size();i++){
+    if (m[i]->nombre().i>0) push_back(m[i]);
+  }
+}
 void Membre::compte(Compteur & c)const{
   for(int i =0; i < size(); i++){
     operator [] (i)->compte(c);
 void Membre::printweight(std::ostream & o) const{
   for(int i=0; i < size(); i++){
     o << operator[](i)->weight();
     if (i < size()-1) o << " ";
   }
+}
 void Membre::coeff( fraction f){
   for (int i=0; i<size(); i++) operator[](i)->coeff(f);
+}
   return result;
+}
 bool Membre::redox()const{
   for (int i=0; i<size(); i++){
     if ((*this)[i]->iselectron()) /* c'est un électron */ return 1;
+  }
   return 0;
+}
 fraction  Membre::nbelectron()const{
   return result;
+}
 void Membre::printNernst(std::ostream & o){
   bool printed = 0;
-  char * prefix="";
+  const char * prefix="";
   for (int i=0; i<size(); i++) {
     if (i>0) prefix="\\,";
     if (operator[](i)->printNernst(o, prefix)){
       printed = 1;
+    }
   return o;
+}
 Membre operator & (Membre & m1, Membre & m2){
   Membre result;
+  result.printnorm(std::cout);
   fraction min(1);
   for(Membre::iterator i = m1.begin(); i < m1.end(); i++){
     for(Membre::iterator j = m2.begin(); j < m2.end(); j++){
       if (**i == **j){
         Molec *m = new Molec(**i);
         diff=(*i)->nb - (*j)->nb;
         m->nb=diff;
+      }
+    }
     result.push_back(m);
   }
   return result;
+}
+bool Chemeq::valdefined()const{
+  return val > MINVAL;
+}
+void Chemeq::addChemeq(const Chemeq * c){
+  if (valdefined() && c->valdefined()){
+    long double e1=enthalpy(), e2=c->enthalpy();
+    fraction n1=nbelectron(), n2=c->nbelectron();
+    long double e = e1+e2;
+    fraction n=n1+n2;
+    if (n.i==0) val=expl(-e/R/T0);
+    else val=-e*n.d/n.i/Faraday;
+  } else {
+    val=MINVAL;
+  }
+  gauche->addMembre(c->gauche);
+  droit->addMembre(c->droit);
+  simplifie(true);
+}
+void Chemeq::subChemeq(const Chemeq * c){
+  if (valdefined() && c->valdefined()){
+    long double e1=enthalpy(), e2=c->enthalpy();
+    fraction n1=nbelectron(), n2=c->nbelectron();
+    long double e = e1-e2;
+    fraction n=n1-n2;
+    if (n.i==0) {
+      val=expl(-e/R/T0);
+    } else{
+      val=-e*n.d/n.i/Faraday;
+    }
+  } else {
+    val=MINVAL;
+  }
+  gauche->addMembre(c->droit);
+  droit->addMembre(c->gauche);
+  simplifie(true);
+}
+long double Chemeq::enthalpy() const{
+  fraction n=nbelectron();
+  if (redox()){
+    return -val*n.i/n.d*Faraday;
+  } else {
+    return -R*T0*logl(val);
+  }
+}
+void Chemeq::simplifie(bool tri=false){
+  Membre communs(*gauche & *droit);
+  if (communs.size()>0){
+    Membre * g, *d;
+    g= new Membre(*gauche - communs);
+    d= new Membre(*droit  - communs);
+    delete gauche;
+    delete droit;
+    gauche=g;
+    droit =d;
+  }
+  gauche->eraseNull();
+  droit->eraseNull();
+  if (tri){
+    numerote();
+    triage();
+  }
 }
 void Chemeq::printnorm(std::ostream & o){
   gauche->printnorm(o);
   o << " -> ";
   droit->printnorm(o);
-  if (val>=0){
+  if (val>MINVAL){
     o << " (";
     if (cste!=std::string("")) o << cste << " = ";
     o << val;
     if (redox()) o << " V";
     o << ")";
 void Chemeq::coeff1(){
   fraction mult = gauche->operator[](0)->nombre();
   mult.inverse();
   gauche->coeff(mult);
   droit->coeff(mult);
-  Membre communs(*gauche & *droit);
-  if (communs.size()>0){
+  simplifie();
-    Membre * g, *d;
+  if (!redox()){
-    g= new Membre(*gauche - communs);
-    d= new Membre(*droit  - communs);
-    delete gauche;
+    val = val*mult.i/mult.d;
-    delete droit;
-    gauche=g;
-    droit =d;
+  }
+}
+void Chemeq::multiply(int num, int den){
+  fraction mult(num,den);
+  gauche->coeff(mult);
+  droit->coeff(mult);
+  simplifie();
   if (!redox()){
     val = val*mult.i/mult.d;
+  }
+}
       o << "}";
+    }
     else {
       droit->printNernst(o);
+    }
-    if (val >=0) {
+    if (val > MINVAL) {
       o << "\\,=\\,";
       if (cste!=std::string("")) o << cste << "\\,=\\,";
-      o << val;
+      o << valeur_latex();
+    }
     else{
       o << "\\,=\\,K";
+    }
+  }
   else{ /* c'est une réaction redox */
     o << "E\\,=\\,";
-    if (val >=0) {
+    if (val > MINVAL) {
       o << val;
+    }
     else{
       o << "E_{0}";
     }
     o << "\\,+\\,\\frac{R\\,T}{";
     o << gauche->nbelectron()+droit->nbelectron() << "\\,F}";
     o << "\\log";
     if (gauche->redox()){ /* c'est une réduction */
       ga=gauche; dr=droit;
+    }
     else{ /* c'est une oxydation */
       ga=droit; dr=gauche;
+    }
         o << "(";
         ga->printNernst(o);
         o << ")";
+      }
+  }
+}
+std::string Chemeq::valeur_latex()const{
+  std::ostringstream so;
+  so << val;
+  std::string s(so.str());
+  std::string::size_type epos=s.find('e',0);
+  if (epos!=std::string::npos){
+    s.erase(epos,1);
+    s.insert(epos,"\\times 10^{");
+    s=s+"}";
+  }
+  return (std::string) s;
+}
 std::ostream & operator << (std::ostream & o, const Chemeq & c){
   o << *c.membregauche() << "\\,\\rightarrow\\," << *c.membredroit();
-  if (c.valeur() >=0) {
+  if (c.valeur() > MINVAL) {
     o << "\\,(";
     if (c.constante()!=std::string("")) o << c.constante() << "\\,=\\,";
-    o << c.valeur();
+    o << c.valeur_latex();
     if (c.redox()) o << " V";
     o << ")";
+  }
   return o;
+}
 fraction operator - (fraction f, fraction g){
   fraction result = fraction(f.i*g.d-g.i*f.d, f.d*g.d);
   result.simplifie();
   return result;
+}
+const fraction & minFraction(const fraction& f1, const fraction &f2){
+  if (f1.i*f2.d > f2.i*f1.d) return f1;
+  else return f2;
+}
 void fraction::simplifie(){
   int maxprem = 23;
   int premiers[]={2,3,5,7,11,13,17,19,23,29};

Subversion Repositories wimsdev

(root)/trunk/wims/src/Misc/chemeq/src/chemeq.cc – Rev 10 → 104