[Bf-blender-cvs] SVN commit: /data/svn/bf-blender [30357] branches/soc-2010-rohith291991: Halfedge DS seems to be working, as well as parametrization.
Rohith B V
rohith291991 at gmail.com
Thu Jul 15 02:59:00 CEST 2010
Revision: 30357
http://projects.blender.org/plugins/scmsvn/viewcvs.php?view=rev&root=bf-blender&revision=30357
Author: rohith291991
Date: 2010-07-15 02:59:00 +0200 (Thu, 15 Jul 2010)
Log Message:
-----------
Halfedge DS seems to be working, as well as parametrization. Further testing/checking has to be done, as well as some minor additions.
Modified Paths:
--------------
branches/soc-2010-rohith291991/intern/comiso/extern/CMesh.h
branches/soc-2010-rohith291991/intern/comiso/extern/uv.h
branches/soc-2010-rohith291991/intern/comiso/intern/uv.cpp
branches/soc-2010-rohith291991/source/blender/modifiers/intern/MOD_quadrangulate.c
Modified: branches/soc-2010-rohith291991/intern/comiso/extern/CMesh.h
===================================================================
--- branches/soc-2010-rohith291991/intern/comiso/extern/CMesh.h 2010-07-15 00:55:31 UTC (rev 30356)
+++ branches/soc-2010-rohith291991/intern/comiso/extern/CMesh.h 2010-07-15 00:59:00 UTC (rev 30357)
@@ -108,6 +108,7 @@
double theta; //Angle that the principle direction field makes with the edge represented by v1-v2
double ut[3]; //Direction of u field
double vt[3]; //Direction of v field
+ int v[3]; //Redundancy
}HFace;
@@ -122,6 +123,10 @@
double uflux;//Flux of u through edge from face[0] to face[1]
double vflux;//Flux of v through edge from face[0] to face[1]
+ double je; //rotational extra 1
+ double ke; //rotational extra 2
+ double u; //Parametrization variable 1 (associated with vertex)
+ double v; //Parametrization variable 2 (associated with vertex)
//Optimization variables
int usign;
int vsign;
Modified: branches/soc-2010-rohith291991/intern/comiso/extern/uv.h
===================================================================
--- branches/soc-2010-rohith291991/intern/comiso/extern/uv.h 2010-07-15 00:55:31 UTC (rev 30356)
+++ branches/soc-2010-rohith291991/intern/comiso/extern/uv.h 2010-07-15 00:59:00 UTC (rev 30357)
@@ -7,7 +7,7 @@
extern "C" {
#endif
-void mi_quad_param(CMesh *mesh,double n_quads);
+void mi_quad_param(HMesh *hm);
void copyHC(HMesh *hm, CMesh *cm);
#ifdef __cplusplus
Modified: branches/soc-2010-rohith291991/intern/comiso/intern/uv.cpp
===================================================================
--- branches/soc-2010-rohith291991/intern/comiso/intern/uv.cpp 2010-07-15 00:55:31 UTC (rev 30356)
+++ branches/soc-2010-rohith291991/intern/comiso/intern/uv.cpp 2010-07-15 00:59:00 UTC (rev 30357)
@@ -7,9 +7,6 @@
#include "Geometry"
#include <math.h>
-
-//DONE
-
void compute_rotated_uv_vars(//rotates variable indices
int var, int rot, //u is 2*var, v is 2*var+1
int &uvar, int &vvar, //rotated variable indices
@@ -21,211 +18,172 @@
vvar = 2*var + (rot+1)%2 ;
}
+void compute_edge_flux(HMesh *hm, bool curl_correction = true) { //with optional curl correction
+
+ //Reset the fluxes
+ int totedge,totface,totvert;
+ totvert=hm->numVerts;
+ totedge=hm->numEdges;
+ totface=hm->numFaces;
+
+ for(int i=0;i<totedge;i++)
+ {
+ hm->edges[i].uflux = 0 ;
+ hm->edges[i].uflux = 0 ;
+ }
+
+ //Compute facet_K2 and fluxes per half edge
+ for(int l=0;l<totface;l++) {
+ //recover the edge vectors
+ HEdge* hi = hm->faces[l].edge ;
+
+ Vector3d p[3] ;
+ for(int i=0; i<3; ++i) {
+ p[i][0] = hi->vertex->co[0] ;
+ p[i][1] = hi->vertex->co[1] ;
+ p[i][2] = hi->vertex->co[2] ;
+ hi = hi->next ;
+ }
+
+ Vector3d e[3] ;
+ for(int i=0; i<3; ++i)
+ e[i] = p[i]-p[(i+2)%3] ;
-/*
-
-void compute_edge_flux(HMesh *mesh, bool curl_correction = false) { //with optional curl correction
-
- int a;
-CVect solution;
-Comiso solver;
-int totedge,totface,totvert;
-totvert=mesh->numVerts;
-totedge=mesh->numEdges;
-totface=mesh->numFaces;
-
-solver.m=totface+totedge;
-solver.n=totface+totedge+1;
-solver.mnnz=(solver.n*8);
-
- //Input attributes
- //MapFacetAttribute<Vector3d> facet_K1(map,"K") ;// u scalar field
- //MapHalfedgeAttribute<int> rot(map,"ls") ;//the rotation transitions
-
- //Output attributes
- //MapFacetAttribute<Vector3d> facet_K2(map,"facet_K2") ;//v scalar field
- //MapHalfedgeAttribute<double> uflux(map,"uflux") ;//flux of u through half edge
- //MapHalfedgeAttribute<double> vflux(map,"vflux") ;//flux of v through half edge
-
- //Optimization variable indices and signs
- //MapHalfedgeAttribute<int> uvar(map) ;
- //MapHalfedgeAttribute<int> vvar(map) ;
- //MapHalfedgeAttribute<int> usign(map) ;
- //MapHalfedgeAttribute<int> vsign(map) ;
-
- //Reset the fluxes
- for(int i=0;i<totedge;i++)
- {
- mesh->edges[i].uflux[0] = 0 ;
- mesh->edges[i].uflux[1] = 0 ;
- mesh->edges[i].vflux[0] = 0 ;
- mesh->edges[i].vflux[1] = 0 ;
- }
-
- //Compute v per facet and fluxes per half edge
- for(int i=0;i<totface;i++) {
- //recover the edge vectors
-
- Vector3d p[3] ;
- for(int j=0; j<3; ++j)
- {
- int v=mesh->faces[i].v[j];
- for(int k=0;k<3;k++)
- {
- p[j][k] = mesh->verts[v].co[k];
- }
-
- }
-
- Vector3d e[3] ;
- for(int j=0; j<3; ++j)
- e[j] = p[j]-p[(j+2)%3] ;
-
//compute u for facet
Vector3d normal = e[1].cross(e[0]);
normal.normalize() ;
Vector3d ut ;
Vector3d x;
- x[0]=mesh->verts[mesh->faces[i].v[1]].co[0]-mesh->verts[mesh->faces[i].v[0]].co[0];
- x[1]=mesh->verts[mesh->faces[i].v[1]].co[1]-mesh->verts[mesh->faces[i].v[0]].co[1];
- x[2]=mesh->verts[mesh->faces[i].v[1]].co[2]-mesh->verts[mesh->faces[i].v[0]].co[2];
+ ut[0]=hm->verts[hm->faces[l].v[1]].co[0]-hm->verts[hm->faces[l].v[0]].co[0];
+ ut[1]=hm->verts[hm->faces[l].v[1]].co[1]-hm->verts[hm->faces[l].v[0]].co[1];
+ ut[2]=hm->verts[hm->faces[l].v[1]].co[2]-hm->verts[hm->faces[l].v[0]].co[2];
- AngleAxis<double> aa(mesh->faces[i].theta, normal);
+ AngleAxis<double> aa(hm->faces[l].theta, normal);
ut=Vector3d(aa*ut);
ut.normalize();
+ hm->faces[l].ut[0]=ut[0];
+ hm->faces[l].ut[1]=ut[1];
+ hm->faces[l].ut[2]=ut[2];
- mesh->faces[i].ut[0]=ut[0];
- mesh->faces[i].ut[1]=ut[1];
- mesh->faces[i].ut[2]=ut[2];
-
Vector3d vt = normal.cross(ut);
vt.normalize();
- mesh->faces[i].vt[0]=vt[0];
- mesh->faces[i].vt[1]=vt[1];
- mesh->faces[i].vt[2]=vt[2];
-
- //compute the flux through each edge
- //hi = fi->halfedge() ; //probably useless
- for(int j=0; j<3; ++j) {
- int t=mesh->faces[i].edges[j];
- mesh->edges[t].uflux[0]=ut.dot(e[j]);
- mesh->edges[t].uflux[1]=ut.dot(e[j]);
- mesh->edges[t].vflux[0]=vt.dot(e[j]);
- mesh->edges[t].vflux[1]=vt.dot(e[j]);
+ hm->faces[l].vt[0]=vt[0];
+ hm->faces[l].vt[1]=vt[1];
+ hm->faces[l].vt[2]=vt[2];
+
+ //compute the flux through each edge
+ hi = hm->faces[l].edge;//probably useless
+ for(int i=0; i<3; ++i) {
+ hi->uflux = ut.dot(e[i]) ;
+ hi->vflux = vt.dot(e[i]) ;
+ hi = hi->next ;
+ }
+ }
+
+ if(!curl_correction) return ;
+ //Curl correction tries to ensure that the u fluxes on the edges of a facet
+ //sum to 0. Intuitively, all the entering flux should also get out.
+
+ //TODO Mainly untested and probably buggy
+
+ //mean of the fluxes per edge
+ /* int nvar = 0 ;
+ FOR_EACH_HALFEDGE(Map,map,hi) {
+ if(Utilz::is_resp(hi) && !hi->is_border_edge()) {
+ //the fluxes wrt the facet of the responsible half edge
+ double uf_in = uflux[hi] ;
+ double vf_in = vflux[hi] ;
+ //the fluxes wrt the opposite facet corrected with rotation
+ int uvar_out, vvar_out, usign_out, vsign_out ;
+ int rotation = (4-rot[hi]+2)%4 ; //+2 since at least the opposite edge
+ compute_rotated_uv_vars(nvar,rotation,uvar_out,vvar_out,usign_out,vsign_out) ;
+ //checks if u and v are switched (rotation of pi/2 or 3pi/2)
+ bool swap = (uvar_out%2)==1 ;
+ double uf_out = swap ? vflux[hi->opposite()] : uflux[hi->opposite()] ;
+ double vf_out = swap ? uflux[hi->opposite()] : vflux[hi->opposite()] ;
+ uf_out *= usign_out ;
+ vf_out *= vsign_out ;
+ //means
+ uflux[hi] = 0.5*(uf_in+uf_out) ;
+ vflux[hi] = 0.5*(vf_in+vf_out) ;
+ //assign opposite
+ uflux[hi->opposite()] = swap ? usign_out*vflux[hi] : usign_out*uflux[hi];
+ vflux[hi->opposite()] = swap ? vsign_out*uflux[hi] : vsign_out*vflux[hi];
+
+ //Assign variable indices
+ uvar[hi] = 2*nvar ;
+ vvar[hi] = 2*nvar+1 ;
+ usign[hi] = 1 ;
+ vsign[hi] = 1 ;
+ uvar[hi->opposite()] = uvar_out ;
+ vvar[hi->opposite()] = vvar_out ;
+ usign[hi->opposite()] = usign_out ;
+ vsign[hi->opposite()] = vsign_out ;
+ ++nvar ;
+ } else if(hi->opposite()->is_border()) {
+ //no opposite value to average, just provide variable index
+ uvar[hi] = 2*nvar ;
+ vvar[hi] = 2*nvar+1 ;
+ usign[hi] = 1 ;
+ vsign[hi] = 1 ;
+ ++nvar ;
+ }
+ }
+
+ //Setup least square system
+ gmm::col_matrix< gmm::wsvector< double > > A ;
+ gmm::resize(A,2*nvar,2*nvar) ;
+ std::vector<double> b(2*nvar,0) ;
+
+ FOR_EACH_HALFEDGE(Map,map,hi) {
+ //These equations state that the resulting flux has to match the initial one
+ if(!hi->is_border()) {
+ if( Utilz::is_resp(hi) || hi->is_border_edge()) {
+ A(uvar[hi],uvar[hi]) += 1 ;
+ A(vvar[hi],vvar[hi]) += 1 ;
+
+ //no need to use the using, the flux is already rotated if necessary
+ b[uvar[hi]] += uflux[hi] ;
+ b[vvar[hi]] += vflux[hi] ;
+ }
+ }
+ }
+ FOR_EACH_FACET(Map,map,fi) {
+ //These equations state that the resulting flux should be 0 on a closed curve
+ Map::Halfedge* hi = fi->halfedge() ;
+ Map::Halfedge* hj = fi->halfedge() ;
+ do {
+ do {
+ //u equations
+ A(uvar[hi],uvar[hj]) += usign[hi]*usign[hj] ;
+ //v equations
+ A(vvar[hi],vvar[hj]) += vsign[hi]*vsign[hj] ;
+ } while(hj != fi->halfedge()) ;
+ } while(hi != fi->halfedge()) ;
+ //The rhs is 0 for these equations
+ }
+
+ //The matrix has to be compressed for this solver
+ gmm::csc_matrix<double> Acsc ;
+ gmm::copy(A,Acsc) ;
+
+ COMISO::MISolver solver;
+ std::vector<int> ids_to_round ;
+ std::vector<double> x(2*nvar,0) ;
+ solver.solve(Acsc, x, b, ids_to_round);
+
+ //set the fluxes to the optimized version
+ FOR_EACH_HALFEDGE(Map,map,hi) {
+ if(!hi->is_border()) {
+ uflux[hi] = usign[hi]*x[uvar[hi]] ;
+ vflux[hi] = vsign[hi]*x[vvar[hi]] ;
+ }
+ }*/
+}
-
- }
- }
-
- if(!curl_correction) return ;
- //Curl correction tries to ensure that the u fluxes on the edges of a facet
- //sum to 0. Intuitively, all the entering flux should also get out.
-
- //TODO Mainly untested and probably buggy
-
- //mean of the fluxes per edge
- int nvar = 0 ;
-
- for(int i=0;i<totedge;i++)
- {
- // TODO later if(Utilz::is_resp(hi) && !hi->is_border_edge())
-
- {
- //the fluxes wrt the facet of the responsible half edge
- double uf_in = uflux[hi] ;
- double vf_in = vflux[hi] ;
- //the fluxes wrt the opposite facet corrected with rotation
- int uvar_out, vvar_out, usign_out, vsign_out ;
- int rotation = (4-rot[hi]+2)%4 ; //+2 since at least the opposite edge
@@ Diff output truncated at 10240 characters. @@
More information about the Bf-blender-cvs
mailing list