Commit 7e251fb0 authored by Simon Boyé's avatar Simon Boyé

Vitelotte: More experiments with FV elements. Weird behaviors.

parent beab52f1
This diff is collapsed.
#!/usr/bin/python3
##################
import sys
from sys import argv
from math import cos, sin, pi
##################
def lerp(alpha, p0, p1):
return [ (1-alpha)*c0 + alpha*c1 for c0, c1 in zip(p0, p1) ]
##################
if len(argv) != 3:
print("Usage:", argv[0], "N_RINGS FILENAME")
exit(1)
filename = argv[2]
nrings = int(argv[1])
out = open(filename, 'w')
sys.stdout = out
print("mvg 1.0")
print("attributes none")
print("v 0 0")
for r in range(1, nrings+1):
for tri in range(6):
alpha = [ tri * pi/3, (tri+1) * pi/3 ]
p = [ [ cos(alpha[i]) * r, sin(alpha[i]) * r ] for i in range(2) ]
for v in range(r):
coords = lerp(v / r, p[0], p[1])
print("v", *coords)
inner = 0
outer = 1
for r in range(1, nrings+1):
first_inner = inner
first_outer = outer
for tri in range(6):
for v in range(r-1):
next_inner = inner + 1 if tri != 5 or v != r-2 else first_inner
print("f", inner, outer, outer+1)
print("f", inner, outer+1, next_inner)
inner += 1
outer += 1
if tri == 5:
print("f", first_inner, outer, first_outer)
else:
print("f", inner, outer, outer + 1)
outer += 1
if r == 1:
inner += 1
# print("f", inner+v, outer+v, outer+v+1)
# print("f", inner+v, outer+v+1, inner+v+1)
# print("f", inner+r-1, outer+r-1, outer+r)
......@@ -148,8 +148,10 @@ void VGMeshWithCurves::setNodesFromCurves(unsigned flags_)
if(isGradientConstraint(pc))
{
Eigen::Vector2f v = (position(toVertex(*hit)) - position(fromVertex(*hit))).normalized();
v = Eigen::Vector2f(-v(1), v(0));
Node gn = addNode(grad * v * (halfedgeOrientation(*hit)? -1: 1));
v = Eigen::Vector2f(-v(1), v(0)) * (halfedgeOrientation(*hit)? -1: 1);
Node gn = addNode(grad * v);
// std::cout << "Point gradient constraint: "
// << v.transpose() << ": " << nodeValue(gn).transpose() << "\n";
halfedgeNode(*hit, EDGE_GRADIENT) = gn;
halfedgeOppositeNode(*hit, EDGE_GRADIENT) = gn;
}
......
......@@ -217,6 +217,8 @@ void Document::solve()
{
std::cout << "Solver error: " << m_fvSolver.errorString() << "\n";
}
exportPlot("plot.obj");
}
emit meshUpdated();
......@@ -397,6 +399,45 @@ void Document::openSaveFinalMeshDialog()
}
void Document::exportPlot(const std::string& filename)
{
std::ofstream out(filename.c_str());
for(Mesh::VertexIterator vit = m_solvedMesh.verticesBegin();
vit != m_solvedMesh.verticesEnd(); ++vit)
{
Mesh::Halfedge h = m_solvedMesh.halfedge(*vit);
Mesh::Halfedge hend = h;
Mesh::Node n;
do
{
n = m_solvedMesh.fromVertexValueNode(h);
h = m_solvedMesh.oppositeHalfedge(m_solvedMesh.prevHalfedge(h));
} while(!m_solvedMesh.isValid(n) && h != hend);
out << "v " << m_solvedMesh.position(*vit).transpose() << " ";
if(m_solvedMesh.isValid(n))
out << m_solvedMesh.nodeValue(n)(0);
else
out << "0";
out << "\n";
}
for(Mesh::FaceIterator fit = m_solvedMesh.facesBegin();
fit != m_solvedMesh.facesEnd(); ++fit)
{
Mesh::Halfedge h = m_solvedMesh.halfedge(*fit);
Mesh::Vertex v0 = m_solvedMesh.toVertex(h);
h = m_solvedMesh.nextHalfedge(h);
Mesh::Vertex v1 = m_solvedMesh.toVertex(h);
h = m_solvedMesh.nextHalfedge(h);
Mesh::Vertex v2 = m_solvedMesh.toVertex(h);
if(m_solvedMesh.isValid(v0) && m_solvedMesh.isValid(v1) && m_solvedMesh.isValid(v2))
out << "f " << v0.idx()+1 << " " << v1.idx()+1 << " " << v2.idx()+1 << "\n";
}
}
SetNodeValue::SetNodeValue(Document *doc, Node node, const NodeValue& value,
bool allowMerge)
: m_document(doc), m_node(node), m_newValue(value),
......
......@@ -113,6 +113,7 @@ public slots:
void openSaveSourceMeshDialog();
void openSaveFinalMeshDialog();
void exportPlot(const std::string& filename);
signals:
void selectionChanged();
......
......@@ -5,17 +5,24 @@ set(CMAKE_INCLUDE_CURRENT_DIR ON)
################################################################################
#############################################
set(surface_mesh_SRCS
surface_mesh.cpp
set(fv_SRCS
fv.cpp
)
add_executable(surface_mesh ${surface_mesh_SRCS})
add_executable(fv ${fv_SRCS})
#############################################
#set(surface_mesh_SRCS
# surface_mesh.cpp
#)
#add_executable(surface_mesh ${surface_mesh_SRCS})
#add_dependencies(buildtests deta_orthogonal_derivatives)
#add_test(deta_orthogonal_derivatives deta_orthogonal_derivatives)
#############################################
set(fem_mesh_SRCS
fem_mesh.cpp
)
add_executable(fem_mesh ${fem_mesh_SRCS})
#set(fem_mesh_SRCS
# fem_mesh.cpp
#)
#add_executable(fem_mesh ${fem_mesh_SRCS})
#add_dependencies(buildtests deta_orthogonal_derivatives)
#add_test(deta_orthogonal_derivatives deta_orthogonal_derivatives)
This diff is collapsed.
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment