Skip to content

Instantly share code, notes, and snippets.

@tehrengruber

tehrengruber/Untitled.ipynb

Created May 21, 2017 12:37
Show Gist options
  • Save tehrengruber/28af7fcdddbb032af1897cff45f3581c to your computer and use it in GitHub Desktop.
Save tehrengruber/28af7fcdddbb032af1897cff45f3581c to your computer and use it in GitHub Desktop.
Download ZIP
Display the source blob
Display the rendered blob
Raw
Untitled.ipynb
{
"cells": [
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"Vec3 Triangle\n",
"Vec2 Edge\n",
"Vec Vertex\n",
"\n",
"cell_types(::Gmsh) = \n",
"\n",
"module Gmsh\n",
"\n",
"@production function read_section(channel::AbstractChannel, ::Section{Nodes}, f::IOStream)::AbstractArray{Vertex, 1}\n",
" number_of_nodes = parse(Int, readline(f))\n",
" produce(:message, :number_vertices, number_of_nodes)\n",
" channel=Channel{Tuple{Int, Vec{3, Float64}}}(100) \n",
" task = @async begin\n",
" info(\"Gmsh reading $(number_of_nodes) vertices asynchronously\")\n",
" for n = 1:number_of_nodes\n",
" line = readline(f)\n",
" s = split(line)\n",
" i = parse(Int, s[1])\n",
" @assert length(s) == 4\n",
" put!(channel (i, Vec{3, Float64}(\n",
" parse(Float64, s[2],\n",
" parse(Float64, s[3],\n",
" parse(Float64, s[4]\n",
" )))\n",
" end\n",
" info(\"Gmsh reading vertices finished\")\n",
" end\n",
" (task, channel)\n",
"end\n",
"\n",
"@production function read_section(channel::AbstractChannel, ::Section{Elements}, f::IOStream)\n",
" n = parse(Int, readline(f))\n",
" produce(:message, :number_elements, number_of_nodes)\n",
" task = @async begin\n",
" info(\"Gmsh reading $(n) elements asynchronously\")\n",
" for n = 1:n\n",
" s = map((x) -> parse(Int, x), split(readline(f)))\n",
" @assert length(s) > 4\n",
" el_number, el_type, number_of_tags = s[1:3]\n",
" el_tags = (\n",
" number_of_tags >= 1 ? s[4] : 0,\n",
" number_of_tags >= 2 ? s[5] : 0,\n",
" number_of_tags >= 3 ? s[6] : 0\n",
" )\n",
" node_number_list = s[4+number_of_tags:end]\n",
" put!(channel, el_number, parse_element(el_type, node_number_list, el_tags))\n",
" elements[el_number] = Element(el_type, node_number_list, el_tags)\n",
" elements_type_count[el_type] += 1\n",
" end\n",
" info(\"Gmsh reading vertices finished\")\n",
" end\n",
" (task, channel)\n",
"end\n",
"\n",
"abstract MeshFunction\n",
"\n",
"type MeshFunction\n",
" indices\n",
" values\n",
"end\n",
"\n",
"function MeshFunction()\n",
" \n",
"end\n",
"\n",
"function transpose{T <: MeshFunction}(mesh_function::T)\n",
" \n",
"end\n",
"\n",
"function init_connectity(::ElementType\"3-node line\", ::ElementType\"2-node line\", triangle_idx, v1, v2, v3)\n",
" # create edges\n",
" put!(mesh.connectivity, →(Dimension{2}, Dimension{1}), triangle_idx, (v1, v2))\n",
" put!(mesh.connectivity, →(Dimension{2}, Dimension{1}), triangle_idx, (v2, v3))\n",
" put!(mesh.connectivity, →(Dimension{2}, Dimension{1}), triangle_idx, (v1, v3))\n",
"end\n",
"\n",
"function init_connectivity(::ElementType\"2-node line\", ::ElementType\"1-node point\", el_idx, v1, v2)\n",
" put!(mesh.connectivity, →(Dimension{1}, Dimension{0}), el_idx, v1)\n",
" put!(mesh.connectivity, →(Dimension{1}, Dimension{0}), el_idx, v2)\n",
"end\n",
"\n",
"function transpose_connectivity()\n",
" \n",
"end\n",
"\n",
"function eltype(itr::→())\n",
"\n",
"function parse_cell(::ElementType\"2-node line\", i1, i2)\n",
" put!(mesh.connectivity, \n",
" →(Dimension{1}, Dimension{0}),\n",
" HalfEdge(i1, # index from\n",
" i2, # index to\n",
" # face to the left\n",
" get(::ElementType\"3-node triangle\", :from, i1, :to, i2))\n",
" )\n",
" put!(mesh.connectivity, →(Dimension{0}, Dimension{1}), vertex_index1, edge_idx)\n",
" put!(mesh.connectivity, →(Dimension{0}, Dimension{1}), vertex_index2, edge_idx)\n",
" HalfEdge(i1, i2)\n",
" HalfEdge()\n",
"end\n",
"\n",
"function parse_cell(::ElementType\"3-node triangle\", triangle_idx, v1, v2, v3)\n",
" \n",
" \n",
" Triangle(v1, v2, v3)\n",
"end"
]
},
{
"cell_type": "code",
"execution_count": 28,
"metadata": {
"collapsed": false
},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"WARNING: Method definition isnull(Any) in module Compat at /home/tehrengruber/.julia/v0.5/Compat/src/Compat.jl:1678 overwritten in module Compat at /home/tehrengruber/.julia/v0.5/Compat/src/Compat.jl:1678.\n",
"WARNING: Method definition redirect_stdin(Function, Any) in module Compat at /home/tehrengruber/.julia/v0.5/Compat/src/Compat.jl:1600 overwritten in module Compat at /home/tehrengruber/.julia/v0.5/Compat/src/Compat.jl:1600.\n",
"WARNING: Method definition redirect_stderr(Function, Any) in module Compat at /home/tehrengruber/.julia/v0.5/Compat/src/Compat.jl:1600 overwritten in module Compat at /home/tehrengruber/.julia/v0.5/Compat/src/Compat.jl:1600.\n",
"WARNING: Method definition take!(Main.Base.AbstractIOBuffer) in module Compat at /home/tehrengruber/.julia/v0.5/Compat/src/Compat.jl:1698 overwritten in module Compat at /home/tehrengruber/.julia/v0.5/Compat/src/Compat.jl:1698.\n",
"WARNING: Method definition redirect_stdout(Function, Any) in module Compat at /home/tehrengruber/.julia/v0.5/Compat/src/Compat.jl:1600 overwritten in module Compat at /home/tehrengruber/.julia/v0.5/Compat/src/Compat.jl:1600.\n",
"WARNING: Method definition convert(Type{Tuple{Vararg{T<:Any, #N<:Any}}}, Real) in module FixedSizeArrays at /home/tehrengruber/.julia/v0.5/FixedSizeArrays/src/conversion.jl:18 overwritten in module FixedSizeArrays at /home/tehrengruber/.julia/v0.5/FixedSizeArrays/src/conversion.jl:18.\n",
"WARNING: static parameter N does not occur in signature for VecGen at /home/tehrengruber/.julia/v0.5/PolytopalKomplex/src/polytopes.jl:99.\n",
"The method will not be callable.\n",
"WARNING: Method definition count() in module Iterators at deprecated.jl:49 overwritten in module Iterators at deprecated.jl:49.\n",
"WARNING: Method definition count(Number) in module Iterators at deprecated.jl:49 overwritten in module Iterators at deprecated.jl:49.\n",
"WARNING: Method definition count(Number, Number) in module Iterators at deprecated.jl:49 overwritten in module Iterators at deprecated.jl:49.\n"
]
}
],
"source": [
"workspace()\n",
"include(\"../gmsh.jl\")\n",
"using PolytopalKomplex\n",
"using Gmsh\n",
"#collect(cells(gmsh, ConvexPolytope\"3-node triangle\"()))\n",
"#display(gmsh)"
]
},
{
"cell_type": "code",
"execution_count": 29,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/plain": [
"Gmsh.GmshData{2,2}((\"2.2\",\"0\",\"8\"),[0.0 1.0 … 0.488281 0.492187; 0.0 0.0 … 0.503906 0.5],Gmsh.Cell[Gmsh.Cell(15,[1],(0,1,0)),Gmsh.Cell(15,[2],(0,2,0)),Gmsh.Cell(15,[3],(0,3,0)),Gmsh.Cell(1,[1,52],(0,1,0)),Gmsh.Cell(1,[52,68],(0,1,0)),Gmsh.Cell(1,[68,35],(0,1,0)),Gmsh.Cell(1,[35,115],(0,1,0)),Gmsh.Cell(1,[115,34],(0,1,0)),Gmsh.Cell(1,[34,99],(0,1,0)),Gmsh.Cell(1,[99,33],(0,1,0)) … Gmsh.Cell(2,[16640,292,16641],(0,5,0)),Gmsh.Cell(2,[16639,16641,2528],(0,5,0)),Gmsh.Cell(2,[4481,16634,16640],(0,5,0)),Gmsh.Cell(2,[16634,161,16640],(0,5,0)),Gmsh.Cell(2,[16634,372,161],(0,5,0)),Gmsh.Cell(2,[16640,161,292],(0,5,0)),Gmsh.Cell(2,[2528,16641,8640],(0,5,0)),Gmsh.Cell(2,[16641,162,8640],(0,5,0)),Gmsh.Cell(2,[16641,292,162],(0,5,0)),Gmsh.Cell(2,[8640,162,355],(0,5,0))],[512,32768,0,0,0,0,0,0,0,0,0,0,0,0,16641,0,0,0,0])"
]
},
"execution_count": 29,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"gmsh = Gmsh.read(\"../meshes/triangle_refined.msh\", mesh_dim=2, world_dim=2)\n"
]
},
{
"cell_type": "code",
"execution_count": 30,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/plain": [
"33283-element Array{Gmsh.Cell,1}:\n",
" Gmsh.Cell(15,[1],(0,1,0)) \n",
" Gmsh.Cell(15,[2],(0,2,0)) \n",
" Gmsh.Cell(15,[3],(0,3,0)) \n",
" Gmsh.Cell(1,[1,52],(0,1,0)) \n",
" Gmsh.Cell(1,[52,68],(0,1,0)) \n",
" Gmsh.Cell(1,[68,35],(0,1,0)) \n",
" Gmsh.Cell(1,[35,115],(0,1,0)) \n",
" Gmsh.Cell(1,[115,34],(0,1,0)) \n",
" Gmsh.Cell(1,[34,99],(0,1,0)) \n",
" Gmsh.Cell(1,[99,33],(0,1,0)) \n",
" Gmsh.Cell(1,[33,107],(0,1,0)) \n",
" Gmsh.Cell(1,[107,32],(0,1,0)) \n",
" Gmsh.Cell(1,[32,98],(0,1,0)) \n",
" ⋮ \n",
" Gmsh.Cell(2,[4481,16640,16639],(0,5,0)) \n",
" Gmsh.Cell(2,[16640,16641,16639],(0,5,0))\n",
" Gmsh.Cell(2,[16640,292,16641],(0,5,0)) \n",
" Gmsh.Cell(2,[16639,16641,2528],(0,5,0)) \n",
" Gmsh.Cell(2,[4481,16634,16640],(0,5,0)) \n",
" Gmsh.Cell(2,[16634,161,16640],(0,5,0)) \n",
" Gmsh.Cell(2,[16634,372,161],(0,5,0)) \n",
" Gmsh.Cell(2,[16640,161,292],(0,5,0)) \n",
" Gmsh.Cell(2,[2528,16641,8640],(0,5,0)) \n",
" Gmsh.Cell(2,[16641,162,8640],(0,5,0)) \n",
" Gmsh.Cell(2,[16641,292,162],(0,5,0)) \n",
" Gmsh.Cell(2,[8640,162,355],(0,5,0)) "
]
},
"execution_count": 30,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"gmsh.cells"
]
},
{
"cell_type": "code",
"execution_count": 31,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"using PolytopalKomplex"
]
},
{
"cell_type": "code",
"execution_count": 32,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/plain": [
"Gmsh.GmshData{2,2}((\"2.2\",\"0\",\"8\"),[0.0 1.0 … 0.488281 0.492187; 0.0 0.0 … 0.503906 0.5],Gmsh.Cell[Gmsh.Cell(15,[1],(0,1,0)),Gmsh.Cell(15,[2],(0,2,0)),Gmsh.Cell(15,[3],(0,3,0)),Gmsh.Cell(1,[1,52],(0,1,0)),Gmsh.Cell(1,[52,68],(0,1,0)),Gmsh.Cell(1,[68,35],(0,1,0)),Gmsh.Cell(1,[35,115],(0,1,0)),Gmsh.Cell(1,[115,34],(0,1,0)),Gmsh.Cell(1,[34,99],(0,1,0)),Gmsh.Cell(1,[99,33],(0,1,0)) … Gmsh.Cell(2,[16640,292,16641],(0,5,0)),Gmsh.Cell(2,[16639,16641,2528],(0,5,0)),Gmsh.Cell(2,[4481,16634,16640],(0,5,0)),Gmsh.Cell(2,[16634,161,16640],(0,5,0)),Gmsh.Cell(2,[16634,372,161],(0,5,0)),Gmsh.Cell(2,[16640,161,292],(0,5,0)),Gmsh.Cell(2,[2528,16641,8640],(0,5,0)),Gmsh.Cell(2,[16641,162,8640],(0,5,0)),Gmsh.Cell(2,[16641,292,162],(0,5,0)),Gmsh.Cell(2,[8640,162,355],(0,5,0))],[512,32768,0,0,0,0,0,0,0,0,0,0,0,0,16641,0,0,0,0])"
]
},
"execution_count": 32,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"gmsh"
]
},
{
"cell_type": "code",
"execution_count": 33,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/plain": [
"32768"
]
},
"execution_count": 33,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"Gmsh.number_of_cells(gmsh, Dim{2})"
]
},
{
"cell_type": "code",
"execution_count": 34,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/plain": [
"1-element Array{Tuple{ConvexPolytope\"3-node triangle\",Int64},1}:\n",
" (ConvexPolytope\"3-node triangle\"(),32768)"
]
},
"execution_count": 34,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"# print some information about the elements\n",
"gmsh\n",
"\n",
"# get all element types of Codim{0}\n",
"cell_types = ((cell_type, c) for (cell_type, c) in number_of_cells_per_type(gmsh) if c>0 && dim(cell_type)==Dim{mesh_dim(gmsh)}())\n",
"collect(cell_types)\n",
" #eltype(cell_types)\n",
"\n",
" #idx_mapping = MeshFunction{MeshIndex{cell_type}}(\n",
" # (idx for (idx, el) in gmsh.elements if el.element_type == cell_type])\n",
" # 1:element_type_count[cell_type]\n",
" #)\n",
" #(for idx, el in gmsh.elements if dim(el) == mesh_dim(gmsh))\n",
" #end\n",
" \n",
" #idx_mapping = MeshFunction{dim(cell_type)}(\n",
" # Array{Vec{subentity_count(::cell_type), VertexIndex}, 1}(),\n",
" # 1:element_type_count[cell_type])\n",
" #for \n",
"#end\n",
"#collect(cell_types)\n",
" \n",
"# end\n",
" # 1. create index mapping\n",
" # 2. push into connectivity array →(Codim{0}, Dim{0})\n",
" \n",
"#el_t_c = ((Gmsh.ElementType{i}(), c) for (i, c) in enumerate(gmsh.element_type_count))\n",
"#codim0_el_types = (el_type for (el_type, count) in el_t_c if (count > 0 && codim(gmsh, el_type) == 0))\n",
"\n",
"# collect(codim0_el_types)\n",
" \n",
"# create a mapping from the original indices to a new set of indices, \n",
"# which is consecutive for elements of equal type\n",
"#map((x, y) -> (x, y), enumerate(gmsh.element_type_count))\n",
"#filter(() -> dim())\n",
"#index_mapping = MeshFunction[]\n",
"#gmsh.element_type_count"
]
},
{
"cell_type": "code",
"execution_count": 35,
"metadata": {
"collapsed": false
},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"WARNING: Method definition eltype(Main.Base.Generator) in module Main at In[8]:2 overwritten in module Main at In[35]:2.\n"
]
},
{
"data": {
"text/plain": [
"eltype (generic function with 89 methods)"
]
},
"execution_count": 35,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"import Base.eltype\n",
"eltype(G::Base.Generator) = Base.code_typed(G.f,(eltype(G.iter),))[1].rettype"
]
},
{
"cell_type": "code",
"execution_count": 75,
"metadata": {
"collapsed": false
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Variables:\n",
" #self#::Base.#collect\n",
" itr::Base.Generator{Filter{Gmsh.##14#16{2},Enumerate{Array{Gmsh.Cell,1}}},Gmsh.##13#15}\n",
" isz::Base.SizeUnknown\n",
" et::Type{Tuple{Int64,Gmsh.Cell}}\n",
" v1::UNION{}\n",
" st::UNION{}\n",
" #temp#::UNION{}\n",
"\n",
"Body:\n",
" begin \n",
" NewvarNode(:(v1::UNION{}))\n",
" NewvarNode(:(st::UNION{}))\n",
" NewvarNode(:(#temp#::UNION{}))\n",
" isz::Base.SizeUnknown = $(Expr(:new, :(Base.SizeUnknown))) # line 299:\n",
" et::Type{Tuple{Int64,Gmsh.Cell}} = $(QuoteNode(Tuple{Int64,Gmsh.Cell})) # line 300:\n",
" (Base.isa)(isz::Base.SizeUnknown,Base.SizeUnknown)::Bool # line 301:\n",
" return $(Expr(:invoke, LambdaInfo for grow_to!(::Array{Tuple{Int64,Gmsh.Cell},1}, ::Base.Generator{Filter{Gmsh.##14#16{2},Enumerate{Array{Gmsh.Cell,1}}},Gmsh.##13#15}), :(Base.grow_to!), :((Core.ccall)(:jl_alloc_array_1d,(Core.apply_type)(Core.Array,Tuple{Int64,Gmsh.Cell},1)::Type{Array{Tuple{Int64,Gmsh.Cell},1}},(Core.svec)(Core.Any,Core.Int)::SimpleVector,Array{Tuple{Int64,Gmsh.Cell},1},0,0,0)::Array{Tuple{Int64,Gmsh.Cell},1}), :(itr))) # line 303: # line 304: # line 305: # line 307: # line 308:\n",
" end::UNION{ARRAY{TUPLE{INT64,GMSH.CELL},1},ARRAY{UNION{},1}}\n"
]
}
],
"source": [
"@code_warntype collect(cells(gmsh, ConvexPolytope\"3-node triangle\"()));"
]
},
{
"cell_type": "code",
"execution_count": 36,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"using FixedSizeArrays"
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {
"collapsed": false
},
"outputs": [
{
"ename": "LoadError",
"evalue": "LoadError: UndefVarError: @ConvexPolytope_str not defined\nwhile loading In[1], in expression starting on line 3",
"output_type": "error",
"traceback": [
"LoadError: UndefVarError: @ConvexPolytope_str not defined\nwhile loading In[1], in expression starting on line 3",
""
]
}
],
"source": [
"mf = nothing\n",
"mf_edges = nothing\n",
"function b{CP}(cp::CP, count::Int)\n",
" mf = MeshFunction(CellIndexRange{typeof(cp)}(1, count), \n",
" Array{CellIndex{ConvexPolytope\"1-node point\"}, 2}(vertex_count(cp), count))\n",
" #cells = ()\n",
" \n",
" #idx_generator = let vc=vertex_count(cp)\n",
" # (idx+vc for idx in mf.indices)\n",
" #end\n",
" #for i in idx_gen\n",
" #mf.values[]\n",
" # println(tu)\n",
" #end\n",
" # add Codim{0} → Dim{0} connectivity\n",
" i=1\n",
" for (ci::Int, cell::Gmsh.Cell) in cells(gmsh, cp)\n",
" #@assert i+2<=length(mf.values)\n",
" @inbounds let v=vertices(cell)\n",
" mf.values[i] = CellIndex{ConvexPolytope\"1-node point\"}(v[1])\n",
" mf.values[i+1] = CellIndex{ConvexPolytope\"1-node point\"}(v[2])\n",
" mf.values[i+2] = CellIndex{ConvexPolytope\"1-node point\"}(v[3])\n",
" end\n",
" i+=3\n",
" end\n",
" \n",
" # read in geometry\n",
" #Geometry(cp, PolytopalComplex{Dim{1}, Dim{2}, calc_type(gmsh)}(())\n",
" \n",
" # allocate oriented faces\n",
" mf_edges = MeshFunction(CellIndexRange{ConvexPolytope\"2-node line\"}(1, count*3), \n",
" Array{CellIndex{ConvexPolytope\"1-node point\"}, 2}(2, count*3))\n",
" \n",
" e = (unsafe_wrap(Array, \n",
" convert(Ptr{Vec{3, Int}}, pointer(mf.values)), \n",
" size(mf.values, 2)))\n",
" \n",
" println(\"sizeof: \", sizeof(e))\n",
" println(e[1:5])\n",
" \n",
" i=1\n",
" for j in 1:size(mf.values, 2)\n",
" let v=Vec{3, CellIndex{ConvexPolytope\"1-node point\"}}(mf.values[1, j], mf.values[2, j], mf.values[3, j])\n",
" #mf_edges[i] = CellVertexList(ConvexPolytope\"3-node triangle\"(), mf[j])\n",
" mf_edges.values[1, i] = v[1]\n",
" mf_edges.values[2, i] = v[2]\n",
" mf_edges.values[1, i+1] = v[2]\n",
" mf_edges.values[2, i+1] = v[3]\n",
" mf_edges.values[1, i+2] = v[3]\n",
" mf_edges.values[2, i+2] = v[1]\n",
" \n",
" #d[(convert(Int, v[1]), convert(Int, v[2]))] = i\n",
" #d[(convert(Int, v[2]), convert(Int, v[3]))] = i+1\n",
" #d[(convert(Int, v[1]), convert(Int, v[3]))] = i+2\n",
"\n",
" i+=3\n",
" end\n",
" end\n",
" \n",
" perm = sortperm(mf_edges.values[1, :])\n",
" sorted_edge_origins = sort(mf_edges.values[1, :])\n",
" \n",
" cc=(c)->convert(Int, c)\n",
" \n",
" nbe=0\n",
" for i in 1:size(mf_edges.values, 2)\n",
" # get all half-edges with origin mf_edges.values[2, i]\n",
" poos_idx = searchsorted(sorted_edge_origins, mf_edges.values[2, i])\n",
" # if exists find the opoosite half edge\n",
" for poo_idx in poos_idx\n",
" let pohe_idx=perm[poo_idx]\n",
" if mf_edges.values[1, pohe_idx] == mf_edges.values[2, i] && \n",
" mf_edges.values[2, pohe_idx] == mf_edges.values[1, i]\n",
" nbe+=1\n",
" end\n",
" end\n",
" end\n",
" end\n",
" \n",
" println(\"non boundary edges: \", nbe/length(mf_edges.indices))\n",
" \n",
" # get neighbouring triangles\n",
" #edges = subentities(ConvexPolytope\"2-node line\", cell_idx, cell)\n",
" #flip_orientation(msh, facet_idx)\n",
" \n",
" \n",
" \n",
" # iterate over all edges and search for edge of opposite direction\n",
" \n",
" # transpose edge conn\n",
" \n",
" # link opposite oriented face\n",
" \n",
" \n",
" \n",
" println(\"edges: \", mf_edges)\n",
" \n",
" #for idx, v in mf\n",
" # for j=1:vertex_count(cp)\n",
" # mf[idx, j] = vertices(cell, j)\n",
" # end\n",
" #end\n",
" println(mf)\n",
" \n",
" return (mf, mf_edges)\n",
"end\n",
"a=()->for (cp::ConvexPolytope, count::Int) in cell_types\n",
" return b(cp, count)\n",
"end\n",
"#@code_warntype a()\n",
"@time result = a()"
]
},
{
"cell_type": "code",
"execution_count": 39,
"metadata": {
"collapsed": false
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"sizeof: 786432\n",
"5-elementArray{FixedSizeArrays.Vec{3,Int64},1}:\n",
" Vec(1,4482,419)\n",
" Vec(4482,4483,419)\n",
" Vec(4482,1474,4483)\n",
" Vec(419,4483,451)\n",
" Vec(1474,4484,4483)\n",
"\n",
"non boundary edges: "
]
},
{
"name": "stderr",
"output_type": "stream",
"text": [
"WARNING: Method definition b(#CP<:Any, Int64) in module Main at In[38]:4 overwritten at In[39]:4.\n"
]
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"0.9947916666666666\n",
"edges: MeshFunction Ptr{Void} @0x00007f0c758f9e90\n",
" mapping 98304 dim PolytopalKomplex.Dim{1}() cells of type ConvexPolytope\"2-node line\"() \n",
" to (2,98304) Array{PolytopalKomplex.CellIndex{ConvexPolytope\"1-node point\"},2} \n",
"\n",
"MeshFunction Ptr{Void} @0x00007f0c772c6a90\n",
" mapping 32768 dim PolytopalKomplex.Dim{2}() cells of type ConvexPolytope\"3-node triangle\"() \n",
" to (3,32768) Array{PolytopalKomplex.CellIndex{ConvexPolytope\"1-node point\"},2} \n",
"\n",
" 0.600277 seconds (733.05 k allocations: 35.300 MB, 1.22% gc time)\n"
]
},
{
"data": {
"text/plain": [
"(MeshFunction Ptr{Void} @0x00007f0c772c6a90\n",
" mapping 32768 dim PolytopalKomplex.Dim{2}() cells of type ConvexPolytope\"3-node triangle\"() \n",
" to (3,32768) Array{PolytopalKomplex.CellIndex{ConvexPolytope\"1-node point\"},2} \n",
",MeshFunction Ptr{Void} @0x00007f0c758f9e90\n",
" mapping 98304 dim PolytopalKomplex.Dim{1}() cells of type ConvexPolytope\"2-node line\"() \n",
" to (2,98304) Array{PolytopalKomplex.CellIndex{ConvexPolytope\"1-node point\"},2} \n",
")"
]
},
"execution_count": 39,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"mf = nothing\n",
"mf_edges = nothing\n",
"function b{CP}(cp::CP, count::Int)\n",
" mf = MeshFunction(CellIndexRange{typeof(cp)}(1, count), \n",
" Array{Cell.VertexList{CP, 2}(vertex_count(cp), count))\n",
" #cells = ()\n",
" \n",
" #idx_generator = let vc=vertex_count(cp)\n",
" # (idx+vc for idx in mf.indices)\n",
" #end\n",
" #for i in idx_gen\n",
" #mf.values[]\n",
" # println(tu)\n",
" #end\n",
" # add Codim{0} → Dim{0} connectivity\n",
" i=1\n",
" for (ci::Int, cell::Gmsh.Cell) in cells(gmsh, cp)\n",
" #@assert i+2<=length(mf.values)\n",
" @inbounds let v=vertices(cell)\n",
" mf.values[i] = CellIndex{ConvexPolytope\"1-node point\"}(v[1])\n",
" mf.values[i+1] = CellIndex{ConvexPolytope\"1-node point\"}(v[2])\n",
" mf.values[i+2] = CellIndex{ConvexPolytope\"1-node point\"}(v[3])\n",
" end\n",
" i+=3\n",
" end\n",
" \n",
" # read in geometry\n",
" #Geometry(cp, PolytopalComplex{Dim{1}, Dim{2}, calc_type(gmsh)}(())\n",
" \n",
" # allocate oriented faces\n",
" mf_edges = MeshFunction(CellIndexRange{ConvexPolytope\"2-node line\"}(1, count*3), \n",
" Array{CellIndex{ConvexPolytope\"1-node point\"}, 2}(2, count*3))\n",
" \n",
" e = (unsafe_wrap(Array, \n",
" convert(Ptr{Vec{3, Int}}, pointer(mf.values)), \n",
" size(mf.values, 2)))\n",
" \n",
" println(\"sizeof: \", sizeof(e))\n",
" println(e[1:5])\n",
" \n",
" i=1\n",
" for j in 1:size(mf.values, 2)\n",
" let v=Vec{3, CellIndex{ConvexPolytope\"1-node point\"}}(mf.values[1, j], mf.values[2, j], mf.values[3, j])\n",
" #mf_edges[i] = CellVertexList(ConvexPolytope\"3-node triangle\"(), mf[j])\n",
" mf_edges.values[1, i] = v[1]\n",
" mf_edges.values[2, i] = v[2]\n",
" mf_edges.values[1, i+1] = v[2]\n",
" mf_edges.values[2, i+1] = v[3]\n",
" mf_edges.values[1, i+2] = v[3]\n",
" mf_edges.values[2, i+2] = v[1]\n",
" \n",
" #d[(convert(Int, v[1]), convert(Int, v[2]))] = i\n",
" #d[(convert(Int, v[2]), convert(Int, v[3]))] = i+1\n",
" #d[(convert(Int, v[1]), convert(Int, v[3]))] = i+2\n",
"\n",
" i+=3\n",
" end\n",
" end\n",
" \n",
" perm = sortperm(mf_edges.values[1, :])\n",
" sorted_edge_origins = sort(mf_edges.values[1, :])\n",
" \n",
" cc=(c)->convert(Int, c)\n",
" \n",
" nbe=0\n",
" for i in 1:size(mf_edges.values, 2)\n",
" # get all half-edges with origin mf_edges.values[2, i]\n",
" poos_idx = searchsorted(sorted_edge_origins, mf_edges.values[2, i])\n",
" # if exists find the opoosite half edge\n",
" for poo_idx in poos_idx\n",
" let pohe_idx=perm[poo_idx]\n",
" if mf_edges.values[1, pohe_idx] == mf_edges.values[2, i] && \n",
" mf_edges.values[2, pohe_idx] == mf_edges.values[1, i]\n",
" nbe+=1\n",
" end\n",
" end\n",
" end\n",
" end\n",
" \n",
" println(\"non boundary edges: \", nbe/length(mf_edges.indices))\n",
" \n",
" # get neighbouring triangles\n",
" #edges = subentities(ConvexPolytope\"2-node line\", cell_idx, cell)\n",
" #flip_orientation(msh, facet_idx)\n",
" \n",
" \n",
" \n",
" # iterate over all edges and search for edge of opposite direction\n",
" \n",
" # transpose edge conn\n",
" \n",
" # link opposite oriented face\n",
" \n",
" \n",
" \n",
" println(\"edges: \", mf_edges)\n",
" \n",
" #for idx, v in mf\n",
" # for j=1:vertex_count(cp)\n",
" # mf[idx, j] = vertices(cell, j)\n",
" # end\n",
" #end\n",
" println(mf)\n",
" \n",
" return (mf, mf_edges)\n",
"end\n",
"a=()->for (cp::ConvexPolytope, count::Int) in cell_types\n",
" return b(cp, count)\n",
"end\n",
"#@code_warntype a()\n",
"@time result = a()"
]
},
{
"cell_type": "code",
"execution_count": 41,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/plain": [
"MeshFunction Ptr{Void} @0x00007f0c772c6a90\n",
" mapping 32768 dim PolytopalKomplex.Dim{2}() cells of type ConvexPolytope\"3-node triangle\"() \n",
" to (3,32768) Array{PolytopalKomplex.CellIndex{ConvexPolytope\"1-node point\"},2} \n"
]
},
"execution_count": 41,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"result[1][CellIndex{ConvexPolytope\"3-node triangle\"}(1)]"
]
},
{
"cell_type": "code",
"execution_count": 45,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/plain": [
"96.0"
]
},
"execution_count": 45,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"length(result[2].values)*sizeof(eltype(result[2].values))/1024/1024"
]
},
{
"cell_type": "code",
"execution_count": 33,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/plain": [
"(516,Gmsh.Cell(2,[1,4482,419],(0,5,0)))"
]
},
"execution_count": 33,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"first(cells(gmsh, ConvexPolytope\"3-node triangle\"()))"
]
},
{
"cell_type": "code",
"execution_count": 24,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/plain": [
"Base.Zip2{StepRange{Int64,Int64},Enumerate{Array{Gmsh.Cell,1}}}(1:3:98302,Enumerate{Array{Gmsh.Cell,1}}(Gmsh.Cell[Gmsh.Cell(15,[1],(0,1,0)),Gmsh.Cell(15,[2],(0,2,0)),Gmsh.Cell(15,[3],(0,3,0)),Gmsh.Cell(1,[1,52],(0,1,0)),Gmsh.Cell(1,[52,68],(0,1,0)),Gmsh.Cell(1,[68,35],(0,1,0)),Gmsh.Cell(1,[35,115],(0,1,0)),Gmsh.Cell(1,[115,34],(0,1,0)),Gmsh.Cell(1,[34,99],(0,1,0)),Gmsh.Cell(1,[99,33],(0,1,0)) … Gmsh.Cell(2,[16640,292,16641],(0,5,0)),Gmsh.Cell(2,[16639,16641,2528],(0,5,0)),Gmsh.Cell(2,[4481,16634,16640],(0,5,0)),Gmsh.Cell(2,[16634,161,16640],(0,5,0)),Gmsh.Cell(2,[16634,372,161],(0,5,0)),Gmsh.Cell(2,[16640,161,292],(0,5,0)),Gmsh.Cell(2,[2528,16641,8640],(0,5,0)),Gmsh.Cell(2,[16641,162,8640],(0,5,0)),Gmsh.Cell(2,[16641,292,162],(0,5,0)),Gmsh.Cell(2,[8640,162,355],(0,5,0))]))"
]
},
"execution_count": 24,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"zip(1:3:length(mf.values), cells(gmsh))"
]
},
{
"cell_type": "code",
"execution_count": 30,
"metadata": {
"collapsed": false
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"(516,Gmsh.Cell(2,[1,4482,419],(0,5,0)))\n",
"(517,Gmsh.Cell(2,[4482,4483,419],(0,5,0)))\n",
"(518,Gmsh.Cell(2,[4482,1474,4483],(0,5,0)))\n",
"(519,Gmsh.Cell(2,[419,4483,451],(0,5,0)))\n",
"(520,Gmsh.Cell(2,[1474,4484,4483],(0,5,0)))\n",
"(521,Gmsh.Cell(2,[4484,4485,4483],(0,5,0)))\n",
"(522,Gmsh.Cell(2,[4484,1475,4485],(0,5,0)))\n",
"(523,Gmsh.Cell(2,[4483,4485,451],(0,5,0)))\n",
"(524,Gmsh.Cell(2,[1474,4486,4484],(0,5,0)))\n",
"(525,Gmsh.Cell(2,[4486,4487,4484],(0,5,0)))\n",
"(516,Gmsh.Cell(2,[1,4482,419],(0,5,0)))\n",
"(517,Gmsh.Cell(2,[4482,4483,419],(0,5,0)))\n",
"(518,Gmsh.Cell(2,[4482,1474,4483],(0,5,0)))\n",
"(519,Gmsh.Cell(2,[419,4483,451],(0,5,0)))\n",
"(520,Gmsh.Cell(2,[1474,4484,4483],(0,5,0)))\n",
"(521,Gmsh.Cell(2,[4484,4485,4483],(0,5,0)))\n",
"(522,Gmsh.Cell(2,[4484,1475,4485],(0,5,0)))\n",
"(523,Gmsh.Cell(2,[4483,4485,451],(0,5,0)))\n",
"(524,Gmsh.Cell(2,[1474,4486,4484],(0,5,0)))\n",
"(525,Gmsh.Cell(2,[4486,4487,4484],(0,5,0)))\n"
]
}
],
"source": [
"for cell in take(cells(gmsh, ConvexPolytope\"3-node triangle\"()), 10)\n",
" println(cell)\n",
" end\n",
"for cell in take(cells(gmsh, ConvexPolytope\"3-node triangle\"()), 10)\n",
" println(cell)\n",
" end"
]
},
{
"cell_type": "code",
"execution_count": 43,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/plain": [
"2×16641 Array{Float64,2}:\n",
" 0.0 1.0 0.0 0.0 0.0 0.0 … 0.488281 0.488281 0.492187\n",
" 0.0 0.0 1.0 0.257813 0.492188 0.476563 0.496094 0.503906 0.5 "
]
},
"execution_count": 43,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"gmsh.nodes"
]
},
{
"cell_type": "code",
"execution_count": 27,
"metadata": {
"collapsed": false
},
"outputs": [
{
"ename": "LoadError",
"evalue": "LoadError: UndefVarError: cell_types not defined\nwhile loading In[27], in expression starting on line 6",
"output_type": "error",
"traceback": [
"LoadError: UndefVarError: cell_types not defined\nwhile loading In[27], in expression starting on line 6",
"",
" in anonymous at ./<missing>:?"
]
}
],
"source": [
"# entities in the gmsh file format have globally unique indicies. we however want\n",
"# indicies of entities to be unique only in their element type\n",
"idx_mapping = zeros(Int, number_of_cells(gmsh))\n",
"# if false\n",
"# as such we create an index mapping from the original indices to new indices\n",
"for cell_type in cell_types\n",
" # create index mapping from the input mesh to \n",
" # a consecutive numbering in the element type\n",
" current_idx = 1\n",
" for (idx, el) in cells(gmsh, cell_type)\n",
" idx_mapping[idx] = current_idx\n",
" current_idx+=1\n",
" end\n",
"end"
]
},
{
"cell_type": "code",
"execution_count": 57,
"metadata": {
"collapsed": false
},
"outputs": [
{
"ename": "LoadError",
"evalue": "LoadError: syntax: unexpected }\nwhile loading In[57], in expression starting on line 2",
"output_type": "error",
"traceback": [
"LoadError: syntax: unexpected }\nwhile loading In[57], in expression starting on line 2",
""
]
}
],
"source": [
"function geometry(msh::Gmsh, filter...)\n",
" (idx, (Geometry(el_type(cell), PolytopalComplex{dim(cell), world_dim()}) for (idx, cell) in cells(msh, filter...)))\n",
"end\n",
"\n",
"function initialize_connectivity(msh::PolytopalComplex, ::Codim{0}, ::Dim{0})\n",
" # determine types (and thus information about the memory layout) of codim 0 cells\n",
" \n",
" \n",
" # determine types of cells with codim > 0\n",
"end"
]
},
{
"cell_type": "code",
"execution_count": 50,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/plain": [
"10×10 Array{Int64,2}:\n",
" 1 11 21 31 41 51 61 71 81 91\n",
" 2 12 22 32 42 52 62 72 82 92\n",
" 3 13 23 33 43 53 63 73 83 93\n",
" 4 14 24 34 44 54 64 74 84 94\n",
" 5 15 25 35 45 55 65 75 85 95\n",
" 6 16 26 36 46 56 66 76 86 96\n",
" 7 17 27 37 47 57 67 77 87 97\n",
" 8 18 28 38 48 58 68 78 88 98\n",
" 9 19 29 39 49 59 69 79 89 99\n",
" 10 20 30 40 50 60 70 80 90 100"
]
},
"execution_count": 50,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"collect(let n=10; (i+(j-1)*n for i=1:n, j=1:n) end)"
]
},
{
"cell_type": "code",
"execution_count": 53,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/plain": [
"10×10 Array{Tuple{Int64,Int64},2}:\n",
" (1,1) (1,2) (1,3) (1,4) (1,5) … (1,7) (1,8) (1,9) (1,10) \n",
" (2,1) (2,2) (2,3) (2,4) (2,5) (2,7) (2,8) (2,9) (2,10) \n",
" (3,1) (3,2) (3,3) (3,4) (3,5) (3,7) (3,8) (3,9) (3,10) \n",
" (4,1) (4,2) (4,3) (4,4) (4,5) (4,7) (4,8) (4,9) (4,10) \n",
" (5,1) (5,2) (5,3) (5,4) (5,5) (5,7) (5,8) (5,9) (5,10) \n",
" (6,1) (6,2) (6,3) (6,4) (6,5) … (6,7) (6,8) (6,9) (6,10) \n",
" (7,1) (7,2) (7,3) (7,4) (7,5) (7,7) (7,8) (7,9) (7,10) \n",
" (8,1) (8,2) (8,3) (8,4) (8,5) (8,7) (8,8) (8,9) (8,10) \n",
" (9,1) (9,2) (9,3) (9,4) (9,5) (9,7) (9,8) (9,9) (9,10) \n",
" (10,1) (10,2) (10,3) (10,4) (10,5) (10,7) (10,8) (10,9) (10,10)"
]
},
"execution_count": 53,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"collect(let n=10; ((i, j) for i=1:n, j=1:n) end)"
]
},
{
"cell_type": "code",
"execution_count": 33,
"metadata": {
"collapsed": false,
"scrolled": true
},
"outputs": [
{
"data": {
"text/plain": [
"2097152-element Array{Gmsh.Cell,1}:\n",
" Gmsh.Cell(2,[1,265218,3425],(0,5,0)) \n",
" Gmsh.Cell(2,[265218,265219,3425],(0,5,0)) \n",
" Gmsh.Cell(2,[265218,69122,265219],(0,5,0)) \n",
" Gmsh.Cell(2,[3425,265219,3587],(0,5,0)) \n",
" Gmsh.Cell(2,[69122,265220,265219],(0,5,0)) \n",
" Gmsh.Cell(2,[265220,265221,265219],(0,5,0)) \n",
" Gmsh.Cell(2,[265220,69123,265221],(0,5,0)) \n",
" Gmsh.Cell(2,[265219,265221,3587],(0,5,0)) \n",
" Gmsh.Cell(2,[69122,265222,265220],(0,5,0)) \n",
" Gmsh.Cell(2,[265222,265223,265220],(0,5,0)) \n",
" Gmsh.Cell(2,[265222,20226,265223],(0,5,0)) \n",
" Gmsh.Cell(2,[265220,265223,69123],(0,5,0)) \n",
" Gmsh.Cell(2,[3587,265221,3424],(0,5,0)) \n",
" ⋮ \n",
" Gmsh.Cell(2,[265217,1050624,1050623],(0,5,0)) \n",
" Gmsh.Cell(2,[1050624,1050625,1050623],(0,5,0))\n",
" Gmsh.Cell(2,[1050624,2308,1050625],(0,5,0)) \n",
" Gmsh.Cell(2,[1050623,1050625,134912],(0,5,0)) \n",
" Gmsh.Cell(2,[265217,1050618,1050624],(0,5,0)) \n",
" Gmsh.Cell(2,[1050618,1189,1050624],(0,5,0)) \n",
" Gmsh.Cell(2,[1050618,3042,1189],(0,5,0)) \n",
" Gmsh.Cell(2,[1050624,1189,2308],(0,5,0)) \n",
" Gmsh.Cell(2,[134912,1050625,527872],(0,5,0)) \n",
" Gmsh.Cell(2,[1050625,1188,527872],(0,5,0)) \n",
" Gmsh.Cell(2,[1050625,2308,1188],(0,5,0)) \n",
" Gmsh.Cell(2,[527872,1188,2674],(0,5,0)) "
]
},
"execution_count": 33,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"# calculate area\n",
"c = collect(cell for (idx, cell) in cells(gmsh, ConvexPolytope\"3-node triangle\"()))\n"
]
},
{
"cell_type": "code",
"execution_count": 56,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/plain": [
"2×1050625 Array{Float64,2}:\n",
" 0.0 1.0 0.0 0.0 0.0 0.0 … 0.498535 0.498535 0.499023\n",
" 0.0 0.0 1.0 0.344727 0.31543 0.317383 0.499512 0.500488 0.5 "
]
},
"execution_count": 56,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"gmsh.nodes"
]
},
{
"cell_type": "code",
"execution_count": 46,
"metadata": {
"collapsed": false
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
" 5.742772 seconds (67.11 M allocations: 1.313 GB, 11.95% gc time)\n"
]
},
{
"data": {
"text/plain": [
"-1.5987211554602254e-13"
]
},
"execution_count": 46,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"calc_area = (area::Float64) -> begin\n",
" for (idx, cell) in enumerate(c)\n",
" # grab vertices\n",
" #coordinates(idx, mesh)\n",
" @inbounds let idx1 = cell.node_number_list[1]::Int,\n",
" idx2 = cell.node_number_list[2]::Int,\n",
" idx3 = cell.node_number_list[3]::Int\n",
" m=Mat{3, 3, Float64}((1., gmsh.nodes[1, idx1], gmsh.nodes[2, idx1]), \n",
" (1., gmsh.nodes[1, idx2], gmsh.nodes[2, idx2]), \n",
" (1., gmsh.nodes[1, idx3], gmsh.nodes[2, idx3]))\n",
" #println(abs(det(m)))\n",
" #break\n",
" area+=abs(det(m))\n",
" end\n",
" #coordinates = collect(gmsh.nodes[1:2, idx] for idx in cell.node_number_list)\n",
" #println(abs(det(Mat{3, 3, Float64}((1., coordinates[1]...), (1., coordinates[2]...), (1., coordinates[3]...)))))\n",
" #println(coordinates)\n",
" end\n",
" 1-area\n",
"end\n",
"@time calc_area(0.0)\n",
"#@code_warntype a(gmsh.cells[1])"
]
},
{
"cell_type": "code",
"execution_count": 16,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/plain": [
"2097152-element Array{Tuple{Int64,Gmsh.Cell},1}:\n",
" (4100,Gmsh.Cell(2,[1,265218,3425],(0,5,0))) \n",
" (4101,Gmsh.Cell(2,[265218,265219,3425],(0,5,0))) \n",
" (4102,Gmsh.Cell(2,[265218,69122,265219],(0,5,0))) \n",
" (4103,Gmsh.Cell(2,[3425,265219,3587],(0,5,0))) \n",
" (4104,Gmsh.Cell(2,[69122,265220,265219],(0,5,0))) \n",
" (4105,Gmsh.Cell(2,[265220,265221,265219],(0,5,0))) \n",
" (4106,Gmsh.Cell(2,[265220,69123,265221],(0,5,0))) \n",
" (4107,Gmsh.Cell(2,[265219,265221,3587],(0,5,0))) \n",
" (4108,Gmsh.Cell(2,[69122,265222,265220],(0,5,0))) \n",
" (4109,Gmsh.Cell(2,[265222,265223,265220],(0,5,0))) \n",
" (4110,Gmsh.Cell(2,[265222,20226,265223],(0,5,0))) \n",
" (4111,Gmsh.Cell(2,[265220,265223,69123],(0,5,0))) \n",
" (4112,Gmsh.Cell(2,[3587,265221,3424],(0,5,0))) \n",
" ⋮ \n",
" (2101240,Gmsh.Cell(2,[265217,1050624,1050623],(0,5,0))) \n",
" (2101241,Gmsh.Cell(2,[1050624,1050625,1050623],(0,5,0)))\n",
" (2101242,Gmsh.Cell(2,[1050624,2308,1050625],(0,5,0))) \n",
" (2101243,Gmsh.Cell(2,[1050623,1050625,134912],(0,5,0))) \n",
" (2101244,Gmsh.Cell(2,[265217,1050618,1050624],(0,5,0))) \n",
" (2101245,Gmsh.Cell(2,[1050618,1189,1050624],(0,5,0))) \n",
" (2101246,Gmsh.Cell(2,[1050618,3042,1189],(0,5,0))) \n",
" (2101247,Gmsh.Cell(2,[1050624,1189,2308],(0,5,0))) \n",
" (2101248,Gmsh.Cell(2,[134912,1050625,527872],(0,5,0))) \n",
" (2101249,Gmsh.Cell(2,[1050625,1188,527872],(0,5,0))) \n",
" (2101250,Gmsh.Cell(2,[1050625,2308,1188],(0,5,0))) \n",
" (2101251,Gmsh.Cell(2,[527872,1188,2674],(0,5,0))) "
]
},
"execution_count": 16,
"metadata": {},
"output_type": "execute_result"
}
],
"source": []
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"using FixedSizeArrays"
]
},
{
"cell_type": "code",
"execution_count": 38,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/plain": [
"FixedSizeArrays.Mat{2,2,Float64}(\n",
" 1.0 1.0\n",
" 2.0 2.0\n",
")\n"
]
},
"execution_count": 38,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"Mat{2, 2, Float64}((1, coordinates[1]...), (1, coordinates[2]...), (1, coordinates[3]...))"
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/plain": [
"PolytopalComplex.Dim{2}"
]
},
"execution_count": 5,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"using PolytopalComplex\n",
"Dim{3}-Dim{2}\n",
"dim(ConvexPolytope\"3-node triangle\"())\n",
"Dim{2}-Codim{0}"
]
},
{
"cell_type": "code",
"execution_count": 11,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"using PolytopalComplex"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/plain": [
"19-element Array{Tuple{Gmsh.ElementType,Int64},1}:\n",
" (Gmsh.ElementType{1}(),4) \n",
" (Gmsh.ElementType{2}(),2) \n",
" (Gmsh.ElementType{3}(),0) \n",
" (Gmsh.ElementType{4}(),0) \n",
" (Gmsh.ElementType{5}(),0) \n",
" (Gmsh.ElementType{6}(),0) \n",
" (Gmsh.ElementType{7}(),0) \n",
" (Gmsh.ElementType{8}(),0) \n",
" (Gmsh.ElementType{9}(),0) \n",
" (Gmsh.ElementType{10}(),0)\n",
" (Gmsh.ElementType{11}(),0)\n",
" (Gmsh.ElementType{12}(),0)\n",
" (Gmsh.ElementType{13}(),0)\n",
" (Gmsh.ElementType{14}(),0)\n",
" (Gmsh.ElementType{15}(),3)\n",
" (Gmsh.ElementType{16}(),0)\n",
" (Gmsh.ElementType{17}(),0)\n",
" (Gmsh.ElementType{18}(),0)\n",
" (Gmsh.ElementType{19}(),0)"
]
},
"execution_count": 4,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"collect((Gmsh.ElementType{i}(), c) for (i, c) in enumerate(gmsh.element_type_count))"
]
},
{
"cell_type": "code",
"execution_count": 10,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/plain": [
"1-element Array{Gmsh.ElementType{2},1}:\n",
" Gmsh.ElementType{2}()"
]
},
"execution_count": 10,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"el_t_c = ((Gmsh.ElementType{i}(), c) for (i, c) in enumerate(gmsh.element_type_count))\n",
" \n",
"collect(el_t_c)\n",
" collect((el_type for (el_type, count) in el_t_c if (count > 0 && dim(el_type) == 2)))"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"\n",
"# initialize connectivity arrays\n",
"connectivity=\n",
"for eltype in eltypes(msh)\n",
" \n",
"end\n",
"\n",
"# map indices\n",
"for idx, el in elements\n",
" if dim(eltype(el))==mesh_dim(mesh)\n",
" connectivity[→(Codim{0}(), Dim{0}())][idx] = vertices(el)\n",
" end\n",
"end\n",
"\n",
"boundary = MeshFunction(filter(is_boundary, mesh[\"1-node line\"]))\n",
"boundary[3]\n",
"@mesh_function (mesh, is_on_boundary) \"1-node line\""
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"c=Channel{Int}(100000)\n",
"d=false\n",
"s=30*1024^2\n",
"producer=@task begin\n",
" tic()\n",
" for i=1:(s/100000)\n",
" @sync for j=1:100000\n",
" push!(c, i)\n",
" end\n",
" end\n",
" d=true\n",
" t=toc()\n",
" println((sizeof(Int)*s)/t/1024/1024)\n",
"end\n",
"\n",
"consumer=@task begin\n",
" while !d\n",
" take!(c)\n",
" end\n",
"end\n",
"\n",
"schedule(producer)\n",
"schedule(consumer)\n"
]
}
],
"metadata": {
"anaconda-cloud": {},
"kernelspec": {
"display_name": "Julia 0.5.1-pre",
"language": "julia",
"name": "julia-0.5"
},
"language_info": {
"file_extension": ".jl",
"mimetype": "application/julia",
"name": "julia",
"version": "0.5.1"
}
},
"nbformat": 4,
"nbformat_minor": 1
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment