MBDyn-1.7.3
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups
genfm.cc File Reference
#include "mbconfig.h"
#include <cmath>
#include "dataman.h"
#include "genfm.h"
#include "bisec.h"
Include dependency graph for genfm.cc:

Go to the source code of this file.

Functions

ElemReadGenericAerodynamicForce (DataManager *pDM, MBDynParser &HP, unsigned int uLabel)
 
static GenericAerodynamicDataReadGenericAerodynamicData (const std::string &fname, doublereal dScaleAngle, doublereal dScaleLength, bool bAlphaFirst)
 
ElemReadGenericAerodynamicForce (DataManager *pDM, MBDynParser &HP, const DofOwner *pDO, unsigned int uLabel)
 

Variables

static const doublereal dAlphaMax [] = { M_PI/2, M_PI }
 
static const doublereal dBetaMax [] = { M_PI, M_PI/2 }
 

Function Documentation

static GenericAerodynamicData* ReadGenericAerodynamicData ( const std::string &  fname,
doublereal  dScaleAngle,
doublereal  dScaleLength,
bool  bAlphaFirst 
)
static

Definition at line 582 of file genfm.cc.

References GenericAerodynamicData::Alpha, GenericAerodynamicData::bAlphaFirst, GenericAerodynamicData::Beta, buf, c, dAlphaMax, GenericAerodynamicData::Data, dBetaMax, LINE_MAX, MBDYN_EXCEPT_ARGS, GenericAerodynamicData::nAlpha, GenericAerodynamicData::name, and GenericAerodynamicData::nBeta.

Referenced by ReadGenericAerodynamicForce().

584 {
585  std::ifstream in(fname.c_str());
586 
587  if (!in) {
588  silent_cerr("ReadGenericAerodynamicData: "
589  "unable to open file \"" << fname << "\""
590  << std::endl);
592  }
593 
594  char buf[LINE_MAX];
595  int nAlpha, nBeta;
596  int c;
597 
598  /* skip comments */
599  for (c = in.get(); c == '%' || c == '#'; c = in.get()) {
600  /* discard to end of line */
601  in.getline(buf, sizeof(buf));
602  }
603  in.putback(c);
604 
605  /* get the size of the matrices */
606  in >> nAlpha >> nBeta;
607  /* discard to end of line */
608  in.getline(buf, sizeof(buf));
609 
610  if (!in) {
611  silent_cerr("ReadGenericAerodynamicData(" << fname << "): "
612  "unable to read size of data matrix "
613  "from file \"" << fname << "\"" << std::endl);
615  }
616 
617  if (nAlpha <= 0) {
618  silent_cerr("ReadGenericAerodynamicData(" << fname << "): "
619  "invalid number of angles of attack " << nAlpha << " "
620  "from file \"" << fname << "\"" << std::endl);
622  }
623 
624  if (nBeta <= 0) {
625  silent_cerr("ReadGenericAerodynamicData(" << fname << "): "
626  "invalid number of sideslip angles " << nBeta << " "
627  "from file \"" << fname << "\"" << std::endl);
629  }
630 
631  /* skip comments */
632  for (c = in.get(); c == '%' || c == '#'; c = in.get()) {
633  /* discard to end of line */
634  in.getline(buf, sizeof(buf));
635  }
636  in.putback(c);
637 
638  if (!in) {
639  silent_cerr("ReadGenericAerodynamicData(" << fname << "): "
640  "unable to get to data "
641  "in file \"" << fname << "\"" << std::endl);
643  }
644 
646  pData->name = fname;
647  pData->bAlphaFirst = bAlphaFirst;
648  pData->nAlpha = nAlpha;
649  pData->nBeta = nBeta;
650 
651  pData->Alpha.resize(nAlpha);
652  pData->Beta.resize(nBeta);
653  pData->Data.resize(nBeta);
654 
655  doublereal dScaleForce = dScaleLength*dScaleLength;
656  doublereal dScaleMoment = dScaleForce*dScaleLength;
657  doublereal dScale[6] = {
658  dScaleForce, dScaleForce, dScaleForce,
659  dScaleMoment, dScaleMoment, dScaleMoment
660  };
661 
662  /* get the matrices */
663  for (int iBeta = 0; iBeta < nBeta; iBeta++) {
664  pData->Data[iBeta].resize(nAlpha);
665 
666  for (int iAlpha = 0; iAlpha < nAlpha; iAlpha++) {
667  doublereal dCoef;
668 
669  /* read (and check) alpha */
670  in >> dCoef;
671  if (iBeta == 0) {
672  if (iAlpha == 0) {
673  doublereal dErr = dCoef*dScaleAngle + ::dAlphaMax[bAlphaFirst];
674  if (std::abs(dErr) > std::numeric_limits<doublereal>::epsilon()) {
675  silent_cerr("ReadGenericAerodynamicData(" << fname << "): "
676  "warning, alpha[0] != -pi/2 (error=" << 100*dErr/::dAlphaMax[bAlphaFirst] << "%)" << std::endl);
677  }
678  } else if (iAlpha == nAlpha - 1) {
679  doublereal dErr = dCoef*dScaleAngle - ::dAlphaMax[bAlphaFirst];
680  if (std::abs(dErr) > std::numeric_limits<doublereal>::epsilon()) {
681  silent_cerr("ReadGenericAerodynamicData(" << fname << "): "
682  "warning, alpha[" << iAlpha << "] != pi/2 (error=" << 100*dErr/::dAlphaMax[bAlphaFirst] << "%)" << std::endl);
683  }
684  }
685 
686  pData->Alpha[iAlpha] = dCoef;
687 
688  } else if (dCoef != pData->Alpha[iAlpha]) {
689  silent_cerr("ReadGenericAerodynamicData"
690  "(" << fname << "): "
691  "inconsistent data, "
692  "Alpha[" << iAlpha << "]"
693  "=" << dCoef << " "
694  "for Beta[" << iBeta << "]"
695  "=" << pData->Beta[iBeta] << " "
696  "differs from previous, "
697  << pData->Alpha[iAlpha]
698  << std::endl);
699  delete pData;
701  }
702 
703  /* read (and check) beta */
704  in >> dCoef;
705  if (iAlpha == 0) {
706  if (iBeta == 0) {
707  doublereal dErr = dCoef*dScaleAngle + ::dBetaMax[bAlphaFirst];
708  if (std::abs(dErr) > std::numeric_limits<doublereal>::epsilon()) {
709  silent_cerr("ReadGenericAerodynamicData(" << fname << "): "
710  "warning, beta[0] != -pi (error=" << 100*dErr/::dBetaMax[bAlphaFirst] << "%)" << std::endl);
711  }
712  } else if (iBeta == nBeta - 1) {
713  doublereal dErr = dCoef*dScaleAngle - ::dBetaMax[bAlphaFirst];
714  if (std::abs(dErr) > std::numeric_limits<doublereal>::epsilon()) {
715  silent_cerr("ReadGenericAerodynamicData(" << fname << "): "
716  "warning, beta[" << iBeta << "] != pi (error=" << 100*dErr/::dBetaMax[bAlphaFirst] << "%)" << std::endl);
717  }
718  }
719 
720  pData->Beta[iBeta] = dCoef;
721 
722  } else if (dCoef != pData->Beta[iBeta]) {
723  silent_cerr("ReadGenericAerodynamicData"
724  "(" << fname << "): "
725  "inconsistent data, "
726  "Beta[" << iBeta << "]"
727  "=" << dCoef << " "
728  "for Alpha[" << iAlpha << "]"
729  "=" << pData->Alpha[iAlpha] << " "
730  "differs from previous, "
731  << pData->Beta[iBeta]
732  << std::endl);
733  delete pData;
735  }
736 
737  for (int iCoef = 0; iCoef < 6; iCoef++) {
738  in >> dCoef;
739 
740  pData->Data[iBeta][iAlpha].dCoef[iCoef] = dCoef*dScale[iCoef];
741  }
742 
743  /* discard to end of line */
744  if (iAlpha < nAlpha - 1 && iBeta < nBeta - 1) {
745  in.getline(buf, sizeof(buf));
746 
747  if (!in) {
748  silent_cerr("ReadGenericAerodynamicData"
749  "(" << fname << "): "
750  "unable to read data past "
751  "iAlpha=" << iAlpha << ", "
752  "iBeta=" << iBeta << std::endl);
753  delete pData;
755  }
756  }
757  }
758  }
759 
760  /* deg => radian */
761  for (int iAlpha = 0; iAlpha < nAlpha; iAlpha++) {
762  pData->Alpha[iAlpha] *= dScaleAngle;
763 
764  if (iAlpha == 0) {
765  continue;
766  }
767 
768  if ( pData->Alpha[iAlpha] <= pData->Alpha[iAlpha - 1]) {
769  silent_cerr("ReadGenericAerodynamicData"
770  "(" << fname << "): "
771  "strict ordering violated between "
772  "Alpha #" << iAlpha - 1 << " and "
773  "Alpha #" << iAlpha << std::endl);
774  delete pData;
776  }
777  }
778 
779  /* deg => radian */
780  for (int iBeta = 0; iBeta < nBeta; iBeta++) {
781  pData->Beta[iBeta] *= dScaleAngle;
782 
783  if (iBeta == 0) {
784  continue;
785  }
786 
787  if ( pData->Beta[iBeta] <= pData->Beta[iBeta - 1]) {
788  silent_cerr("ReadGenericAerodynamicData"
789  "(" << fname << "): "
790  "strict ordering violated between "
791  "Beta #" << iBeta - 1 << " and "
792  "Beta #" << iBeta << std::endl);
793  delete pData;
795  }
796  }
797 
798  return pData;
799 }
static const doublereal dAlphaMax[]
Definition: genfm.cc:44
#define MBDYN_EXCEPT_ARGS
Definition: except.h:63
std::vector< doublereal > Alpha
Definition: genfm.h:52
std::string name
Definition: genfm.h:44
std::vector< std::vector< GenericAerodynamicCoef > > Data
Definition: genfm.h:66
std::vector< doublereal > Beta
Definition: genfm.h:53
#define LINE_MAX
Definition: mbdyn.h:163
static const doublereal dBetaMax[]
Definition: genfm.cc:45
static std::stack< cleanup * > c
Definition: cleanup.cc:59
static doublereal buf[BUFSIZE]
Definition: discctrl.cc:333
double doublereal
Definition: colamd.c:52
Elem* ReadGenericAerodynamicForce ( DataManager pDM,
MBDynParser HP,
unsigned int  uLabel 
)
Elem* ReadGenericAerodynamicForce ( DataManager pDM,
MBDynParser HP,
const DofOwner pDO,
unsigned int  uLabel 
)

Definition at line 802 of file genfm.cc.

References Elem::AERODYNAMIC, GenericAerodynamicData::bAlphaFirst, Eye3, DataManager::fReadOutput(), IncludeParser::GetFileName(), IncludeParser::GetLineData(), MBDynParser::GetPosRel(), HighParser::GetReal(), MBDynParser::GetRotRel(), HighParser::GetYesNo(), HighParser::IsKeyWord(), M_PI, MBDYN_EXCEPT_ARGS, ReadGenericAerodynamicData(), DataManager::ReadNode(), SAFENEWWITHCONSTRUCTOR, Node::STRUCTURAL, and Zero3.

804 {
805  /* Nodo */
806  const StructNode* pNode = pDM->ReadNode<const StructNode, Node::STRUCTURAL>(HP);
807 
808  /* The offset is in the reference frame of the node */
809  ReferenceFrame RF(pNode);
810  Vec3 f(Zero3);
811  if (HP.IsKeyWord("position")) {
812  f = HP.GetPosRel(RF);
813  }
814 
815  /* the orientation is in flight mechanics (FIXME?) reference frame:
816  * X forward, Y to the right, Z down */
817  Mat3x3 Ra(Eye3);
818  if (HP.IsKeyWord("orientation")) {
819  Ra = HP.GetRotRel(RF);
820  }
821 
822  /* 1. by default, which means that coefficients are only normalized
823  * by the dynamic pressure */
824  doublereal dRefSurface = 1.;
825  doublereal dRefLength = 1.;
826  bool bAlphaFirst(false);
827 
828  if (HP.IsKeyWord("reference" "surface")) {
829  dRefSurface = HP.GetReal();
830  }
831 
832  if (HP.IsKeyWord("reference" "length")) {
833  dRefLength = HP.GetReal();
834  }
835 
836  if (HP.IsKeyWord("alpha" "first")) {
837  if (!HP.GetYesNo(bAlphaFirst)) {
838  silent_cerr("GenericAerodynamicForce(" << uLabel << "): "
839  "invalid \"alpha first\" value at line " << HP.GetLineData() << std::endl);
841  }
842  }
843 
844  /* TODO: allow to reference previously loaded data */
845  GenericAerodynamicData *pData = 0;
846  if (HP.IsKeyWord("file")) {
847  doublereal dScaleAngle = 1.;
848  if (HP.IsKeyWord("angle" "units")) {
849  if (HP.IsKeyWord("radians")) {
850  dScaleAngle = 1.;
851 
852  } else if (HP.IsKeyWord("degrees")) {
853  dScaleAngle = M_PI/180.;
854 
855  } else {
856  silent_cerr("GenericAerodynamicForce(" << uLabel << "): "
857  "unknown angle unit at line " << HP.GetLineData() << std::endl);
859  }
860 
861  } else if (HP.IsKeyWord("scale" "angles")) {
862  doublereal d = HP.GetReal(dScaleAngle);
863  if (d <= 0.) {
864  silent_cerr("GenericAerodynamicForce(" << uLabel << "): "
865  "invalid angle scale factor " << d
866  << " at line " << HP.GetLineData() << std::endl);
868  }
869  dScaleAngle = d;
870  }
871 
872  doublereal dScaleLength = 1.;
873  if (HP.IsKeyWord("scale" "lengths")) {
874  doublereal d = HP.GetReal(dScaleLength);
875  if (d <= 0.) {
876  silent_cerr("GenericAerodynamicForce(" << uLabel << "): "
877  "invalid length scale factor " << d
878  << " at line " << HP.GetLineData() << std::endl);
880  }
881  dScaleLength = d;
882  }
883 
884  std::string fname(HP.GetFileName());
885  pData = ReadGenericAerodynamicData(fname, dScaleAngle, dScaleLength, bAlphaFirst);
886 
887  } else if (HP.IsKeyWord("reference")) {
888  silent_cerr("GenericAerodynamicForce(" << uLabel << "): "
889  "references to generic aerodynamic data not implemented yet "
890  "at line " << HP.GetLineData() << std::endl);
892 
893  } else {
894  silent_cerr("GenericAerodynamicForce(" << uLabel << "): "
895  "keyword \"file\" or \"reference\" expected "
896  "at line " << HP.GetLineData() << std::endl);
898  }
899 
900  flag fOut = pDM->fReadOutput(HP, Elem::AERODYNAMIC);
901 
902  Elem *pEl = 0;
904  GenericAerodynamicForce(uLabel, pDO,
905  pNode, f, Ra,
906  dRefSurface, dRefLength, bAlphaFirst,
907  pData, fOut));
908 
909  return pEl;
910 }
flag fReadOutput(MBDynParser &HP, const T &t) const
Definition: dataman.h:1064
Mat3x3 GetRotRel(const ReferenceFrame &rf)
Definition: mbpar.cc:1795
#define M_PI
Definition: gradienttest.cc:67
const Vec3 Zero3(0., 0., 0.)
long int flag
Definition: mbdyn.h:43
#define MBDYN_EXCEPT_ARGS
Definition: except.h:63
Definition: matvec3.h:98
const Mat3x3 Eye3(1., 0., 0., 0., 1., 0., 0., 0., 1.)
virtual const char * GetFileName(enum Delims Del=DEFAULTDELIM)
Definition: parsinc.cc:673
static GenericAerodynamicData * ReadGenericAerodynamicData(const std::string &fname, doublereal dScaleAngle, doublereal dScaleLength, bool bAlphaFirst)
Definition: genfm.cc:582
Vec3 GetPosRel(const ReferenceFrame &rf)
Definition: mbpar.cc:1331
virtual bool IsKeyWord(const char *sKeyWord)
Definition: parser.cc:910
#define SAFENEWWITHCONSTRUCTOR(pnt, item, constructor)
Definition: mynewmem.h:698
virtual bool GetYesNo(bool &bRet)
Definition: parser.cc:1022
Definition: elem.h:75
double doublereal
Definition: colamd.c:52
virtual HighParser::ErrOut GetLineData(void) const
Definition: parsinc.cc:697
Node * ReadNode(MBDynParser &HP, Node::Type type) const
Definition: dataman3.cc:2309
virtual doublereal GetReal(const doublereal &dDefval=0.0)
Definition: parser.cc:1056

Here is the call graph for this function:

Variable Documentation

const doublereal dAlphaMax[] = { M_PI/2, M_PI }
static

Definition at line 44 of file genfm.cc.

Referenced by GenericAerodynamicForce::AssVec(), and ReadGenericAerodynamicData().

const doublereal dBetaMax[] = { M_PI, M_PI/2 }
static

Definition at line 45 of file genfm.cc.

Referenced by GenericAerodynamicForce::AssVec(), and ReadGenericAerodynamicData().