ffglut.cpp

FreeFem++可以生成高质量的有限元网格。可以用于流体力学

    
    Bmin= M - D2;
    Bmax= M + D2;
    
    if((debug > 10)) cout << " cadreortho: "<< " :: " << Bmin << " " << Bmax <<" oB " << oBmin << " " << oBmax << endl;
    
    // if((debug > 10)) cout << "cadreortho\n";
}
void OneWindow::setLighting()
{
  if(withlight)
    {
      if(plotdim==3)
	{
	  GLfloat lp0[4] = { cam.x,cam.y,cam.z, 1.0 };
	  glLightfv(GL_LIGHT0,GL_POSITION,lp0);	
	  
	  if(debug>1)  cout << " Light pos  3d:  " << cam << endl;
	}
      else
	{
	  
	  GLfloat position[] = {Pvue3.x,Pvue3.y,Pvue3.z+(Bmax3.z-Bmin3.z)*3,1.f} ;
	  glLightfv(GL_LIGHT0, GL_POSITION, position);
	  
	}
      
      float cca=0.3,ccd=1., ccs=0.8;
      GLfloat ambient[] = {cca,cca,cca,1.0f};//diff閞ents param鑤res
      GLfloat diffuse[] = {ccd,ccd,ccd,1.0f};
      GLfloat specular_reflexion[] = {ccs,ccs,ccs,1.0f};
      GLubyte shiny_obj = 128;
      glLightModeli(GL_LIGHT_MODEL_TWO_SIDE,1);
      glEnable(GL_LIGHTING);//positionnement de la lumi鑢e avec
      glLightfv(GL_LIGHT0,GL_AMBIENT,ambient);//les diff閞ents param鑤res
      glLightfv(GL_LIGHT0,GL_DIFFUSE,diffuse);
      
      glEnable(GL_COLOR_MATERIAL);//sp閏ification de la r閒lexion sur les mat閞iaux
      glColorMaterial(GL_FRONT_AND_BACK,GL_AMBIENT_AND_DIFFUSE);
      glMaterialfv(GL_FRONT_AND_BACK,GL_AMBIENT,ambient);
      glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,diffuse);
      // glMaterialfv(GL_FRONT_AND_BACK,GL_SPECULAR,specular_reflexion);// on peut le faire avant chaque objet
      //glMateriali(GL_FRONT_AND_BACK,GL_SHININESS,shiny_obj);//si on veut qu'ils aient des caract閞istiques #
      glShadeModel(GL_FLAT);  
      glEnable(GL_LIGHTING);
      glEnable(GL_LIGHT0); 
    }
  else
    {
      glDisable(GL_LIGHTING);
      glDisable(GL_LIGHT0);
    }

}

void OneWindow::unsetLighting()
{
  glDisable(GL_LIGHTING);
  glDisable(GL_LIGHT0);
  //g->lightning=false;
}

OnePlotBorder::OnePlotBorder(PlotStream & f) 
  :OnePlot(4,2,1)
{ 
    long nbd;
    f>> nbd;
    data.resize(nbd);
    for(int i=0;i<nbd;++i)
      {
	  long n; 
	  f>> n;
	  //cout << n << endl;
	  data[i].resize(n+1);
	  for(int j=0;j<=n;++j)
	    {
		long l; 
		double x,y;
		f >> l>> x >> y;
		// cout << x << ' ' << y << ' ' <<  l << endl;
		R2 P(x,y);
		Pmin=Minc(Pmin,P);
		Pmax=Maxc(Pmax,P);
		data[i][j]=make_pair(l,P);
	    }
      }
    ffassert(f.good());
}

void OnePlot::GLDraw(OneWindow *win)
{
    ThePlot & plot= *win->theplot; 
    //      glNewList(gllists,GL_COMPILE_AND_EXECUTE); // save  la list sans affichage
      Draw(win);
      //  glEndList();  // fin de la list
    //  else glCallList(gllists);
}

void ThePlot::DrawHelp(OneWindow *win) 
{
  int i = 1;
  win->Show("Enter a keyboard character in the FreeFem Graphics window in order to:",i++);
  
  i+=1;
  win->Show("enter) wait next plot",i++);
  win->Show("p)     previous plot",i++);
  win->Show("ESC)   exit from ffglut",i++);
  win->Show("?)  show this help window",i++);
  win->Show("+) -)   zoom in/out  around the cursor 3/2 times ",i++);
  win->Show("=)   reset vue ",i++);
  win->Show("r)   refresh plot ",i++);
  win->Show("3)   switch 3d/2d plot (in test)  keys : ",i++);
  win->Show("        move : <- ->  mouse to return etc.  ",i++);
  win->Show("        z) Z) (focal zoom unzoom)  ",i++);
  win->Show("        H) h) switch increase or decrease the Z scale of the plot ",i++);

  win->Show("a) A) increase or decrease the size arrow ",i++);
  win->Show("b)  switch between black and white or color plotting ",i++);
  win->Show("g)  switch between grey or color plotting ",i++);
  win->Show("f)  switch between filling iso or not  ",i++);
  win->Show("l)  switch between lighting or not  ",i++);
  win->Show("v)  switch between show the numerical value of iso or not",i++);
  win->Show("w)  graphic window dump in file ",i++);
  win->Show("m)  switch between show  meshes or not",i++);
  win->Show("w)  window dump in file ffglutXXXX.ppm ",i++);
  win->Show("any other key : nothing ",++i);
}


void ThePlot::Draw(OneWindow *win) 
{
    if((debug>1 )) 
      { 
	cout << "      ThePlot::Plot " << count << " " << this << " win " << win << " " << state ;
	for (list<OnePlot *>::iterator i= plots.begin();i != plots.end(); ++i)
	  cout << (**i).what;
	cout << endl;
	
      }
    if(state==0) {
      state=1;
      win->DefaultView();
    }

    win->SetView();
    for (list<OnePlot *>::iterator i= plots.begin();i != plots.end(); ++i)
	(*i)->Draw(win);
    
    if(cm || value)
      { //  screen plot ...
	win->SetScreenView();
	if(cm)
	  {
	    color(1);
	    win->DrawCommentaire(cm->c_str(),0.1,0.97);
	  }
	if(value)
	  {  
	    int k0=0;
	    if(withiso)
	      {
		win->PlotValue(Viso,k0,"IsoValue");
		k0=1+Viso.N();
	      }
	    if(witharrow)
	      win->PlotValue(Varrow,k0,"Vec Value");
	  }
	// for picking..
	win->SetView();
      }
    changeViso=false;
    changeVarrow=false;
    changeColor=false;
    changeBorder=false;
    changeFill=false;
}

void  ThePlot::SetColorTable(int nb)
{
    tbc.resize(nb);
    for (int i=0;i<nb;++i)
	tbc[i].set(i,nb,this);    
}



ThePlot::ThePlot(PlotStream & fin,ThePlot *old,int kcount)
  :  count(kcount), state(0),
     changeViso(true),changeVarrow(true),changeColor(true),
     changeBorder(true),changeFill(true), withiso(false),witharrow(false),
     plotdim(2),theta(30.*M_PI/180.),phi(20.*M_PI/180.),dcoef(1),focal(20.*M_PI/180.),
     datadim(1)
     
{
  
  hsv=true; // hsv  type   
  coeff=1;
  wait=0;
  value=false;
  fill=false;
  aspectratio=false;
  clean=true;
  uaspectratio=false;
  pViso=false;
  pVarrow=false;
  Niso=0;
  Narrow=20;
  bw=false;
  psfile=0;
  cm=0;
  grey=0;
  if(old) {
    grey=old->grey;
  }
  greyo=grey;
  drawborder=true;
  drawmeshes=false;
  add=false; 
  keepPV=false;
    
  Pmin=R3(+dinfty,+dinfty,+dinfty);
  fmin = +dinfty;    
  fmax = -dinfty;
  Pmax=R3(-dinfty,-dinfty,-dinfty);
  vmax=0;
  
  coefr=1;
  long dimpp=0;
  long cas;
  
  while(1)
    {
      fin >> cas;
      if((debug > 4)) cout << " read cas: " << cas << "  " << PlotStream::dt_endarg << endl;
      
      if(cas==PlotStream::dt_endarg) break;
      switch (cas) {
      case  0: fin >> coeff; break;
      case 1: fin >> cm; break;
      case 2: fin >> psfile; break;
      case 3: fin >> wait; break;
      case 4: fin >> fill; break;
      case 5: fin >> value; break;
      case 6: fin >> clean; break;
      case 7: fin >> aspectratio;uaspectratio=true; break;
      case 8: fin >> boundingbox; break;
      case 9: fin >> Niso; break;
      case 10: fin >> Narrow; break;
      case 11: fin >> Viso;Niso=Viso.N();pViso=true; break;
      case 12: fin >> Varrow;Narrow=Varrow.N();pVarrow=true; break;
      case 13: fin >> bw; break;
      case 14: fin >> grey; break;
      case 15: fin >> colors; break;
      case 16: fin >> drawborder; break;
      case 17: fin >> dimpp; break;// ajout fevr 2008  v3.0.6
      case 18: fin >> add; break;
      case 19: fin >> keepPV; break;
      default: 
	cout << " cas : " << cas <<endl;
	ffassert(0);
	break;
      }
      ffassert(fin.good());
    }
  if(dimpp) plotdim=dimpp; 
  //    if( !uaspectratio) aspectratio= true;
  ffassert(cas==PlotStream::dt_endarg);
  if((debug > 2))
    {
      cout << " coeff " << coeff <<", ";
      if(cm)
	cout << " cm " << *cm <<", ";
      if(wait)	cout << " wait " ;
      if(fill)	cout << " fill " ;
      if(value)	cout << " value " ;
      if(bw)	cout << " bw " ;
      if(grey)	cout << " grey " ;
      if(drawborder)	cout << " drawborder " ;
      if(colors.N()) cout << "\n colors =" << colors;
      if(boundingbox.N()) cout << "\n bb  =" << boundingbox;
      
      cout << endl;
    } 
  fin.GetMeshes();  
  long nbmeshes;
  fin >> nbmeshes;
  if((debug > 2)) cout << " read nb : mesh " << nbmeshes << endl;
  Ths.resize(nbmeshes);
  for(int i=0;i<nbmeshes;++i)
    Ths[i]=0;
  for(int i=0;i<nbmeshes;++i)
    { 
      long l;
      fin >> l;
      if(l>=0) 
	{
	if((debug > 3)) cout << " read mesh " << i  << " -> " << l << "  " <<nbmeshes << endl;
	l--;
	ffassert(l>=0 && l < nbmeshes);
	ffassert(Ths[l]==0);
	fin >>Ths[l] ;
	if((debug > 3))
	  cout << i << " nt/nv " << l << " "  <<Ths[l]->nt << " " << Ths[l]->nv << endl;
	ffassert(fin.good());
      }
      else // Add FH optimisation FH 11/12/2008 (not use to day)
	{// the mesh is already in the previous plot with number ll
	  ffassert(l==-1);
	  long ll;
	  fin >> l>> ll; // read l and ll
	  ffassert(old);
	  Ths[l]=old->Ths[ll];
	  Ths[l]->increment(); // 
	}
      
    }
    long nbmeshes3=0;
   if(fin.GetMeshes3())
     { //  il y a des solution 3d; 
	 
	 fin >> nbmeshes3;
	 if((debug > 2)) cout << " read nb : mesh3 " << nbmeshes3 << endl;
	 Ths3.resize(nbmeshes3);
	 for(int i=0;i<nbmeshes3;++i)
	     Ths3[i]=0;
	 for(int i=0;i<nbmeshes3;++i)
	   { 
	       long l;
	       fin >> l;
	       if(l>=0) 
		 {
		     if((debug > 3)) cout << " read mesh3 " << i  << " -> " << l << "  " <<nbmeshes3 << endl;
		     l--;
		     ffassert(l>=0 && l < nbmeshes3);
		     ffassert(Ths3[l]==0);
		     fin >>Ths3[l] ;
		     if((debug > 3))
			 cout << i << " nt/nv " << l << " "  <<Ths3[l]->nt << " " << Ths3[l]->nv << endl;
		     ffassert(fin.good());
		 }
	       else // Add FH optimisation FH 11/12/2008 (not use to day)
		 {// the mesh is already in the previous plot with number ll
		     ffassert(l==-1);
		     long ll;
		     fin >> l>> ll; // read l and ll
		     ffassert(old);
		     Ths3[l]=old->Ths3[ll];
		     Ths3[l]->increment(); // 
		 }
	       
	   }	 
	 
	 
	 
      fin.GetPlots();
     }
  long nbplot;
  int iso3d=0;
  fin >>nbplot;
  if((debug > 2)) cout << " nb item plot " << nbplot << endl;
  for(int i=0;i<nbplot;++i)
    {
      long what;
      OnePlot *p=0;
      fin >> what;
      long imsh;
      if((debug > 2)) cout << "    plot  " << i << " what " << what << endl;
      if(what !=3 && !uaspectratio) aspectratio= true;
      if(what==-1)  // gestion of error (empty plot)
	  p = new OnePlotError(fin);
      else if(what==0) 
	{ 
	  
	  fin >> imsh;
	  p=new OnePlotMesh(Ths[imsh-1]);
	}
      else if (what==1 || what==2)
	{
	  fin >> imsh;
	  if(what==1) withiso=true;
	  else if (what==2) witharrow=true;
	  if((debug > 10)) cout << " plot : mesh " << imsh << endl;
	  ffassert(imsh>0 && imsh <=nbmeshes);
	  p=new OnePlotFE(Ths[imsh-1],what,fin);
	}
      else if(what==3)
	p=new OnePlotCurve(fin);
      else if(what==4)
	p=new OnePlotBorder(fin);
      else if(what==5) 
	{ 	    
	    fin >> imsh;
	    p=new OnePlotMesh3(Ths3[imsh-1]);
	}
      else if (what==6 )
	{
	  iso3d++;
	    fin >> imsh;
	    if(what==6) withiso=true;
	   
	    if((debug > 10)) cout << " plot : mesh3 " << imsh << endl;
	    ffassert(imsh>0 && imsh <=nbmeshes3);
	    p=new OnePlotFE3(Ths3[imsh-1],what,fin);
	}
	
      else
	{
	  cout << "Bizarre unkown what :  " << what<< endl;
	ffassert(0);
	}
      ffassert(p);
      plots.push_back(p);
      p->bb(Pmin,Pmax);
      p->bfv(fmin,fmax,vmax);
      plotdim=max(plotdim,p->dim);
      ffassert(fin.good());		      
      datadim=max(datadim,p->dim); 
    }
  if(Niso==0) 
      Niso = iso3d ? 5 : 20;
    
  // cout << "\t\t\t\t  f min, max v max :" << fmin << " " << fmax << " " << vmax << endl;
    
  double ref_f = abs(fmax)+abs(fmin) ; 
  if(fmax < fmin)
    {
      fmax = 1;
      fmin = 0;
    }
  else if( (fmax-fmin) <= 1e-8*ref_f)
    {
      if(ref_f< 1e-20) ref_f=0.5;
      fmax += ref_f/2;
      fmin -= ref_f/2;
    }
  PminT=Pmin;
  PmaxT=Pmax;
  fminT=fmin;
  fmaxT=fmax;
  if(old && 0)
    {
      Pmin= Minc(Pmin,old->PminT);
      Pmax= Maxc(Pmax,old->PmaxT);
      fmax= Max(fmax,old->fmaxT);
      fmin= Min(fmin,old->fminT);
    }
  
  z0= fminT +(fmaxT-fminT)*0.01;
  if((debug > 2)) cout << "               data bound: " << PminT << " " << PmaxT  << " fmin == " << fminT << "  " << fmaxT 
		       << " z0 " << z0 <<  endl;
  fin.GetEndPlot(); 
  Viso.resize(Niso);
  Varrow.resize(Narrow);
  
  SetColorTable(Max(Niso,Narrow)+4) ;           
  
}


void ThePlot::SetDefIsoV()
{
  R d = fill ? (fmaxT-fminT)/(Niso-2)  : (fmaxT-fminT)/(Niso-1);       
  R x = fill ? fminT-d/2 :fminT+d/2;
  if(!pViso)
    {
      for (int i = 0;i < Niso;i++)
	{Viso[i]=x;x +=d; }
      if (fill ) {Viso[0]=fminT-d;Viso[Niso-1]=fmaxT+d;}    
    }
  x=0; 
  if(debug>10)
    cout << "vmax=  " << vmax << endl;
  d= sqrt(vmax)/(Narrow-1.1);   
  if (!pVarrow)
    for (int i = 0;i < Narrow;i++)
      {Varrow[i]=x;x +=d; }
  SetColorTable(Max(Niso,Narrow)+4) ; 
}

void OneWindow::Show(const char *str,int i)
{
  int hx= 18;
  int ix= width/20;
  int iy= height-hx*i;
  plot(ix,iy,str,3);
}

void  FillRectRasterPos(R x0,R y0,R x1,R y1)
{
  //  if((debug > 10)) cout << "FR Rp:   " << x0 << " " << y0 << " " << x1 << " " << y1 << endl;
  glPolygonMode(GL_FRONT_AND_BACK,GL_FILL);//GL_FILL
  glBegin(GL_POLYGON);
  glVertex2d(x0, y0);
  glVertex2d(x1, y0);
  glVertex2d(x1, y1);
  glVertex2d(x0, y1);
  glEnd();
  
}

void  OneWindow::FillRect(R x0,R y0,R x1,R y1)
{
  
  double z1=  (zmin+zmax)/2; // line 
  glPolygonMode(GL_FRONT,GL_FILL);//GL_FILL
  glBegin(GL_POLYGON);
  glVertex3d(x0, y0,z1);
  glVertex3d(x1, y0,z1);
  glVertex3d(x1, y1,z1);
  glVertex3d(x0, y1,z1);
  glEnd();
}

void OneWindow::PlotValue(const KN_<double> & Viso,int  k0,const char * cmm)
{
  
  ShowGlerror("PlotValue b");
  //  glRasterPos2f(x,y);
  if((debug > 10)) cout << "PlotValue:" << cmm << " " << k0 << " " << width << " " <<height << endl;
  R xmin=0,xmax=width,ymin=0,ymax=height;
  if((debug > 10)) cout << "PlotValue " << Viso << endl;
  // int ix,iy;
  // GetSizeScreen(ix,iy);   
  
  R dx=(xmax-xmin);
  R dy=(ymax-ymin);
  //  10 points  
  // int kk = Max(30,iy/10);
  R h=10;
  R ho=h*1.1;
  R x0=xmin+dx*0.85;
  R y= ymin+dy*0.97;
  if((debug > 10)) cout << x0 << " " << y << " " << h <<  endl;
  y -= k0* ho;
  this->color(0);
  FillRectRasterPos(x0-h*0.5,y-h*(1.4*Viso.N()+0.3),x0+h*9,y+h*1.5);
  ShowGlerror("PlotValue m");
  this->color(1);

  plot(x0+ho,y,cmm);
  y -=  ho;
  for (int i=0;i<Viso.N();i++)
    {
      if((debug > 10)) cout << " x0 : " << x0<< " " << y << " " << h << " v = " << Viso[i] << endl;
      this->color(i+4);
      FillRectRasterPos(x0,y,x0+h,y+h);
      plot(x0+ho,y+3*h/10,Viso[i]);
      y -=  ho;
      ;                 
    }
  ShowGlerror("PlotValue f");
}


void OneWindow::DrawCommentaire(const char * cm,R x,R y) 
{
  
    R xmin=0,xmax=height,ymin=0,ymax=height;
    float dx=(xmax-xmin);
    float dy=(ymax-ymin);    
    plot(xmin+dx*x,ymin+dy*y,cm);   
}

/*