dolfin team mailing list archive
-
dolfin team
-
Mailing list archive
-
Message #18017
Re: Merge request for some sandbox stuff
Could you publish your sandbox branch?
I think we should just merge published sandbox branches without (or
with minimal) review, but it's good to use the possibility of
publishing branches on Launchpad + I need to practice how to
review/approve/merge branches on Launchpad. I've never done it before.
--
Anders
On Tue, Apr 06, 2010 at 11:25:07AM +0200, Andre Massing wrote:
> Hi!
>
> While I am working on a implementation of Nitsche methods for overlapping
> meshes in DOLFIN I produce some sand, which might be nice to be added to the
> sandbox :) It contains among other things some prototypes for building a
> complete intersection map for two overlapping meshes (which cells/facets are
> covered by which cells/facest) and for barycenter quadrature formulas on
> arbitrary polyhedrons and polygons (unfortunately just for the 3D case, up to
> now).
>
> Kind regards,
> Andre
> # Bazaar merge directive format 2 (Bazaar 0.90)
> # revision_id: massing@xxxxxxxxx-20100406091417-ogh30d38n95eh5cq
> # target_branch: bzr+ssh://bazaar.launchpad.net/~dolfin-\
> # core/dolfin/sandbox/
> # testament_sha1: 74c2e1fd30a903bacbd8c50606cb683173e73813
> # timestamp: 2010-04-06 11:15:17 +0200
> # base_revision_id: logg@xxxxxxxxx-20100226124738-bk9jlyb9voycqbwa
> #
> # Begin patch
> === renamed directory 'cgal' => 'nitsche_method'
> === added directory 'nitsche_method/data'
> === added directory 'nitsche_method/data/meshes'
> === added file 'nitsche_method/data/meshes/rotator.xml'
> --- nitsche_method/data/meshes/rotator.xml 1970-01-01 00:00:00 +0000
> +++ nitsche_method/data/meshes/rotator.xml 2010-03-16 14:39:09 +0000
> @@ -0,0 +1,402 @@
> +<?xml version="1.0" encoding="UTF-8"?>
> +
> +<dolfin xmlns:dolfin="http://www.fenics.org/dolfin/";>
> + <mesh celltype="triangle" dim="2">
> + <vertices size="117">
> + <vertex index="0" x="-0.1" y="-1"/>
> + <vertex index="1" x="0" y="-1"/>
> + <vertex index="2" x="0" y="-0.9"/>
> + <vertex index="3" x="-0.1" y="-0.9"/>
> + <vertex index="4" x="0.1" y="-1"/>
> + <vertex index="5" x="0.1" y="-0.9"/>
> + <vertex index="6" x="0" y="-0.8"/>
> + <vertex index="7" x="-0.1" y="-0.8"/>
> + <vertex index="8" x="0.1" y="-0.8"/>
> + <vertex index="9" x="0" y="-0.7"/>
> + <vertex index="10" x="-0.1" y="-0.7"/>
> + <vertex index="11" x="0.1" y="-0.7"/>
> + <vertex index="12" x="0" y="-0.6"/>
> + <vertex index="13" x="-0.1" y="-0.6"/>
> + <vertex index="14" x="0.1" y="-0.6"/>
> + <vertex index="15" x="0" y="-0.5"/>
> + <vertex index="16" x="-0.1" y="-0.5"/>
> + <vertex index="17" x="0.1" y="-0.5"/>
> + <vertex index="18" x="0" y="-0.4"/>
> + <vertex index="19" x="-0.1" y="-0.4"/>
> + <vertex index="20" x="0.1" y="-0.4"/>
> + <vertex index="21" x="0" y="-0.3"/>
> + <vertex index="22" x="-0.1" y="-0.3"/>
> + <vertex index="23" x="0.1" y="-0.3"/>
> + <vertex index="24" x="0" y="-0.2"/>
> + <vertex index="25" x="-0.1" y="-0.2"/>
> + <vertex index="26" x="0.1" y="-0.2"/>
> + <vertex index="27" x="0" y="-0.1"/>
> + <vertex index="28" x="-0.1" y="-0.1"/>
> + <vertex index="29" x="0.1" y="-0.1"/>
> + <vertex index="30" x="-1" y="-0.1"/>
> + <vertex index="31" x="-0.9" y="-0.1"/>
> + <vertex index="32" x="-0.9" y="0"/>
> + <vertex index="33" x="-1" y="0"/>
> + <vertex index="34" x="-0.8" y="-0.1"/>
> + <vertex index="35" x="-0.8" y="0"/>
> + <vertex index="36" x="-0.7" y="-0.1"/>
> + <vertex index="37" x="-0.7" y="0"/>
> + <vertex index="38" x="-0.6" y="-0.1"/>
> + <vertex index="39" x="-0.6" y="0"/>
> + <vertex index="40" x="-0.5" y="-0.1"/>
> + <vertex index="41" x="-0.5" y="0"/>
> + <vertex index="42" x="-0.4" y="-0.1"/>
> + <vertex index="43" x="-0.4" y="0"/>
> + <vertex index="44" x="-0.3" y="-0.1"/>
> + <vertex index="45" x="-0.3" y="0"/>
> + <vertex index="46" x="-0.2" y="-0.1"/>
> + <vertex index="47" x="-0.2" y="0"/>
> + <vertex index="48" x="-0.1" y="0"/>
> + <vertex index="49" x="0" y="0"/>
> + <vertex index="50" x="0.1" y="0"/>
> + <vertex index="51" x="0.2" y="-0.1"/>
> + <vertex index="52" x="0.2" y="0"/>
> + <vertex index="53" x="0.3" y="-0.1"/>
> + <vertex index="54" x="0.3" y="0"/>
> + <vertex index="55" x="0.4" y="-0.1"/>
> + <vertex index="56" x="0.4" y="0"/>
> + <vertex index="57" x="0.5" y="-0.1"/>
> + <vertex index="58" x="0.5" y="0"/>
> + <vertex index="59" x="0.6" y="-0.1"/>
> + <vertex index="60" x="0.6" y="0"/>
> + <vertex index="61" x="0.7" y="-0.1"/>
> + <vertex index="62" x="0.7" y="0"/>
> + <vertex index="63" x="0.8" y="-0.1"/>
> + <vertex index="64" x="0.8" y="0"/>
> + <vertex index="65" x="0.9" y="-0.1"/>
> + <vertex index="66" x="0.9" y="0"/>
> + <vertex index="67" x="1" y="-0.1"/>
> + <vertex index="68" x="1" y="0"/>
> + <vertex index="69" x="-0.9" y="0.1"/>
> + <vertex index="70" x="-1" y="0.1"/>
> + <vertex index="71" x="-0.8" y="0.1"/>
> + <vertex index="72" x="-0.7" y="0.1"/>
> + <vertex index="73" x="-0.6" y="0.1"/>
> + <vertex index="74" x="-0.5" y="0.1"/>
> + <vertex index="75" x="-0.4" y="0.1"/>
> + <vertex index="76" x="-0.3" y="0.1"/>
> + <vertex index="77" x="-0.2" y="0.1"/>
> + <vertex index="78" x="-0.1" y="0.1"/>
> + <vertex index="79" x="0" y="0.1"/>
> + <vertex index="80" x="0.1" y="0.1"/>
> + <vertex index="81" x="0.2" y="0.1"/>
> + <vertex index="82" x="0.3" y="0.1"/>
> + <vertex index="83" x="0.4" y="0.1"/>
> + <vertex index="84" x="0.5" y="0.1"/>
> + <vertex index="85" x="0.6" y="0.1"/>
> + <vertex index="86" x="0.7" y="0.1"/>
> + <vertex index="87" x="0.8" y="0.1"/>
> + <vertex index="88" x="0.9" y="0.1"/>
> + <vertex index="89" x="1" y="0.1"/>
> + <vertex index="90" x="0" y="0.2"/>
> + <vertex index="91" x="-0.1" y="0.2"/>
> + <vertex index="92" x="0.1" y="0.2"/>
> + <vertex index="93" x="0" y="0.3"/>
> + <vertex index="94" x="-0.1" y="0.3"/>
> + <vertex index="95" x="0.1" y="0.3"/>
> + <vertex index="96" x="0" y="0.4"/>
> + <vertex index="97" x="-0.1" y="0.4"/>
> + <vertex index="98" x="0.1" y="0.4"/>
> + <vertex index="99" x="0" y="0.5"/>
> + <vertex index="100" x="-0.1" y="0.5"/>
> + <vertex index="101" x="0.1" y="0.5"/>
> + <vertex index="102" x="0" y="0.6"/>
> + <vertex index="103" x="-0.1" y="0.6"/>
> + <vertex index="104" x="0.1" y="0.6"/>
> + <vertex index="105" x="0" y="0.7"/>
> + <vertex index="106" x="-0.1" y="0.7"/>
> + <vertex index="107" x="0.1" y="0.7"/>
> + <vertex index="108" x="0" y="0.8"/>
> + <vertex index="109" x="-0.1" y="0.8"/>
> + <vertex index="110" x="0.1" y="0.8"/>
> + <vertex index="111" x="0" y="0.9"/>
> + <vertex index="112" x="-0.1" y="0.9"/>
> + <vertex index="113" x="0.1" y="0.9"/>
> + <vertex index="114" x="0" y="1"/>
> + <vertex index="115" x="-0.1" y="1"/>
> + <vertex index="116" x="0.1" y="1"/>
> + </vertices>
> + <cells size="152">
> + <triangle index="0" v0="0" v1="1" v2="2"/>
> + <triangle index="1" v0="0" v1="2" v2="3"/>
> + <triangle index="2" v0="1" v1="4" v2="5"/>
> + <triangle index="3" v0="1" v1="2" v2="5"/>
> + <triangle index="4" v0="2" v1="3" v2="6"/>
> + <triangle index="5" v0="3" v1="6" v2="7"/>
> + <triangle index="6" v0="2" v1="5" v2="8"/>
> + <triangle index="7" v0="2" v1="6" v2="8"/>
> + <triangle index="8" v0="6" v1="7" v2="9"/>
> + <triangle index="9" v0="7" v1="9" v2="10"/>
> + <triangle index="10" v0="6" v1="8" v2="11"/>
> + <triangle index="11" v0="6" v1="9" v2="11"/>
> + <triangle index="12" v0="9" v1="10" v2="12"/>
> + <triangle index="13" v0="10" v1="12" v2="13"/>
> + <triangle index="14" v0="9" v1="11" v2="14"/>
> + <triangle index="15" v0="9" v1="12" v2="14"/>
> + <triangle index="16" v0="12" v1="13" v2="15"/>
> + <triangle index="17" v0="13" v1="15" v2="16"/>
> + <triangle index="18" v0="12" v1="14" v2="17"/>
> + <triangle index="19" v0="12" v1="15" v2="17"/>
> + <triangle index="20" v0="15" v1="16" v2="18"/>
> + <triangle index="21" v0="16" v1="18" v2="19"/>
> + <triangle index="22" v0="15" v1="17" v2="20"/>
> + <triangle index="23" v0="15" v1="18" v2="20"/>
> + <triangle index="24" v0="18" v1="19" v2="21"/>
> + <triangle index="25" v0="19" v1="21" v2="22"/>
> + <triangle index="26" v0="18" v1="20" v2="23"/>
> + <triangle index="27" v0="18" v1="21" v2="23"/>
> + <triangle index="28" v0="21" v1="22" v2="24"/>
> + <triangle index="29" v0="22" v1="24" v2="25"/>
> + <triangle index="30" v0="21" v1="23" v2="26"/>
> + <triangle index="31" v0="21" v1="24" v2="26"/>
> + <triangle index="32" v0="24" v1="25" v2="27"/>
> + <triangle index="33" v0="25" v1="27" v2="28"/>
> + <triangle index="34" v0="24" v1="26" v2="29"/>
> + <triangle index="35" v0="24" v1="27" v2="29"/>
> + <triangle index="36" v0="30" v1="31" v2="32"/>
> + <triangle index="37" v0="30" v1="32" v2="33"/>
> + <triangle index="38" v0="31" v1="34" v2="35"/>
> + <triangle index="39" v0="31" v1="32" v2="35"/>
> + <triangle index="40" v0="34" v1="36" v2="37"/>
> + <triangle index="41" v0="34" v1="35" v2="37"/>
> + <triangle index="42" v0="36" v1="38" v2="39"/>
> + <triangle index="43" v0="36" v1="37" v2="39"/>
> + <triangle index="44" v0="38" v1="40" v2="41"/>
> + <triangle index="45" v0="38" v1="39" v2="41"/>
> + <triangle index="46" v0="40" v1="42" v2="43"/>
> + <triangle index="47" v0="40" v1="41" v2="43"/>
> + <triangle index="48" v0="42" v1="44" v2="45"/>
> + <triangle index="49" v0="42" v1="43" v2="45"/>
> + <triangle index="50" v0="44" v1="46" v2="47"/>
> + <triangle index="51" v0="44" v1="45" v2="47"/>
> + <triangle index="52" v0="28" v1="46" v2="48"/>
> + <triangle index="53" v0="46" v1="47" v2="48"/>
> + <triangle index="54" v0="27" v1="28" v2="49"/>
> + <triangle index="55" v0="28" v1="48" v2="49"/>
> + <triangle index="56" v0="27" v1="29" v2="50"/>
> + <triangle index="57" v0="27" v1="49" v2="50"/>
> + <triangle index="58" v0="29" v1="51" v2="52"/>
> + <triangle index="59" v0="29" v1="50" v2="52"/>
> + <triangle index="60" v0="51" v1="53" v2="54"/>
> + <triangle index="61" v0="51" v1="52" v2="54"/>
> + <triangle index="62" v0="53" v1="55" v2="56"/>
> + <triangle index="63" v0="53" v1="54" v2="56"/>
> + <triangle index="64" v0="55" v1="57" v2="58"/>
> + <triangle index="65" v0="55" v1="56" v2="58"/>
> + <triangle index="66" v0="57" v1="59" v2="60"/>
> + <triangle index="67" v0="57" v1="58" v2="60"/>
> + <triangle index="68" v0="59" v1="61" v2="62"/>
> + <triangle index="69" v0="59" v1="60" v2="62"/>
> + <triangle index="70" v0="61" v1="63" v2="64"/>
> + <triangle index="71" v0="61" v1="62" v2="64"/>
> + <triangle index="72" v0="63" v1="65" v2="66"/>
> + <triangle index="73" v0="63" v1="64" v2="66"/>
> + <triangle index="74" v0="65" v1="67" v2="68"/>
> + <triangle index="75" v0="65" v1="66" v2="68"/>
> + <triangle index="76" v0="32" v1="33" v2="69"/>
> + <triangle index="77" v0="33" v1="69" v2="70"/>
> + <triangle index="78" v0="32" v1="35" v2="71"/>
> + <triangle index="79" v0="32" v1="69" v2="71"/>
> + <triangle index="80" v0="35" v1="37" v2="72"/>
> + <triangle index="81" v0="35" v1="71" v2="72"/>
> + <triangle index="82" v0="37" v1="39" v2="73"/>
> + <triangle index="83" v0="37" v1="72" v2="73"/>
> + <triangle index="84" v0="39" v1="41" v2="74"/>
> + <triangle index="85" v0="39" v1="73" v2="74"/>
> + <triangle index="86" v0="41" v1="43" v2="75"/>
> + <triangle index="87" v0="41" v1="74" v2="75"/>
> + <triangle index="88" v0="43" v1="45" v2="76"/>
> + <triangle index="89" v0="43" v1="75" v2="76"/>
> + <triangle index="90" v0="45" v1="47" v2="77"/>
> + <triangle index="91" v0="45" v1="76" v2="77"/>
> + <triangle index="92" v0="47" v1="48" v2="78"/>
> + <triangle index="93" v0="47" v1="77" v2="78"/>
> + <triangle index="94" v0="48" v1="49" v2="79"/>
> + <triangle index="95" v0="48" v1="78" v2="79"/>
> + <triangle index="96" v0="49" v1="50" v2="80"/>
> + <triangle index="97" v0="49" v1="79" v2="80"/>
> + <triangle index="98" v0="50" v1="52" v2="81"/>
> + <triangle index="99" v0="50" v1="80" v2="81"/>
> + <triangle index="100" v0="52" v1="54" v2="82"/>
> + <triangle index="101" v0="52" v1="81" v2="82"/>
> + <triangle index="102" v0="54" v1="56" v2="83"/>
> + <triangle index="103" v0="54" v1="82" v2="83"/>
> + <triangle index="104" v0="56" v1="58" v2="84"/>
> + <triangle index="105" v0="56" v1="83" v2="84"/>
> + <triangle index="106" v0="58" v1="60" v2="85"/>
> + <triangle index="107" v0="58" v1="84" v2="85"/>
> + <triangle index="108" v0="60" v1="62" v2="86"/>
> + <triangle index="109" v0="60" v1="85" v2="86"/>
> + <triangle index="110" v0="62" v1="64" v2="87"/>
> + <triangle index="111" v0="62" v1="86" v2="87"/>
> + <triangle index="112" v0="64" v1="66" v2="88"/>
> + <triangle index="113" v0="64" v1="87" v2="88"/>
> + <triangle index="114" v0="66" v1="68" v2="89"/>
> + <triangle index="115" v0="66" v1="88" v2="89"/>
> + <triangle index="116" v0="78" v1="79" v2="90"/>
> + <triangle index="117" v0="78" v1="90" v2="91"/>
> + <triangle index="118" v0="79" v1="80" v2="92"/>
> + <triangle index="119" v0="79" v1="90" v2="92"/>
> + <triangle index="120" v0="90" v1="91" v2="93"/>
> + <triangle index="121" v0="91" v1="93" v2="94"/>
> + <triangle index="122" v0="90" v1="92" v2="95"/>
> + <triangle index="123" v0="90" v1="93" v2="95"/>
> + <triangle index="124" v0="93" v1="94" v2="96"/>
> + <triangle index="125" v0="94" v1="96" v2="97"/>
> + <triangle index="126" v0="93" v1="95" v2="98"/>
> + <triangle index="127" v0="93" v1="96" v2="98"/>
> + <triangle index="128" v0="96" v1="97" v2="99"/>
> + <triangle index="129" v0="97" v1="99" v2="100"/>
> + <triangle index="130" v0="96" v1="98" v2="101"/>
> + <triangle index="131" v0="96" v1="99" v2="101"/>
> + <triangle index="132" v0="99" v1="100" v2="102"/>
> + <triangle index="133" v0="100" v1="102" v2="103"/>
> + <triangle index="134" v0="99" v1="101" v2="104"/>
> + <triangle index="135" v0="99" v1="102" v2="104"/>
> + <triangle index="136" v0="102" v1="103" v2="105"/>
> + <triangle index="137" v0="103" v1="105" v2="106"/>
> + <triangle index="138" v0="102" v1="104" v2="107"/>
> + <triangle index="139" v0="102" v1="105" v2="107"/>
> + <triangle index="140" v0="105" v1="106" v2="108"/>
> + <triangle index="141" v0="106" v1="108" v2="109"/>
> + <triangle index="142" v0="105" v1="107" v2="110"/>
> + <triangle index="143" v0="105" v1="108" v2="110"/>
> + <triangle index="144" v0="108" v1="109" v2="111"/>
> + <triangle index="145" v0="109" v1="111" v2="112"/>
> + <triangle index="146" v0="108" v1="110" v2="113"/>
> + <triangle index="147" v0="108" v1="111" v2="113"/>
> + <triangle index="148" v0="111" v1="112" v2="114"/>
> + <triangle index="149" v0="112" v1="114" v2="115"/>
> + <triangle index="150" v0="111" v1="113" v2="116"/>
> + <triangle index="151" v0="111" v1="114" v2="116"/>
> + </cells>
> + <data>
> + <data_entry name="global vertex indices">
> + <meshfunction type="uint" dim="0" size="117">
> + <entity index="0" value="429"/>
> + <entity index="1" value="430"/>
> + <entity index="2" value="471"/>
> + <entity index="3" value="470"/>
> + <entity index="4" value="431"/>
> + <entity index="5" value="472"/>
> + <entity index="6" value="512"/>
> + <entity index="7" value="511"/>
> + <entity index="8" value="513"/>
> + <entity index="9" value="553"/>
> + <entity index="10" value="552"/>
> + <entity index="11" value="554"/>
> + <entity index="12" value="594"/>
> + <entity index="13" value="593"/>
> + <entity index="14" value="595"/>
> + <entity index="15" value="635"/>
> + <entity index="16" value="634"/>
> + <entity index="17" value="636"/>
> + <entity index="18" value="676"/>
> + <entity index="19" value="675"/>
> + <entity index="20" value="677"/>
> + <entity index="21" value="717"/>
> + <entity index="22" value="716"/>
> + <entity index="23" value="718"/>
> + <entity index="24" value="758"/>
> + <entity index="25" value="757"/>
> + <entity index="26" value="759"/>
> + <entity index="27" value="799"/>
> + <entity index="28" value="798"/>
> + <entity index="29" value="800"/>
> + <entity index="30" value="789"/>
> + <entity index="31" value="790"/>
> + <entity index="32" value="831"/>
> + <entity index="33" value="830"/>
> + <entity index="34" value="791"/>
> + <entity index="35" value="832"/>
> + <entity index="36" value="792"/>
> + <entity index="37" value="833"/>
> + <entity index="38" value="793"/>
> + <entity index="39" value="834"/>
> + <entity index="40" value="794"/>
> + <entity index="41" value="835"/>
> + <entity index="42" value="795"/>
> + <entity index="43" value="836"/>
> + <entity index="44" value="796"/>
> + <entity index="45" value="837"/>
> + <entity index="46" value="797"/>
> + <entity index="47" value="838"/>
> + <entity index="48" value="839"/>
> + <entity index="49" value="840"/>
> + <entity index="50" value="841"/>
> + <entity index="51" value="801"/>
> + <entity index="52" value="842"/>
> + <entity index="53" value="802"/>
> + <entity index="54" value="843"/>
> + <entity index="55" value="803"/>
> + <entity index="56" value="844"/>
> + <entity index="57" value="804"/>
> + <entity index="58" value="845"/>
> + <entity index="59" value="805"/>
> + <entity index="60" value="846"/>
> + <entity index="61" value="806"/>
> + <entity index="62" value="847"/>
> + <entity index="63" value="807"/>
> + <entity index="64" value="848"/>
> + <entity index="65" value="808"/>
> + <entity index="66" value="849"/>
> + <entity index="67" value="809"/>
> + <entity index="68" value="850"/>
> + <entity index="69" value="872"/>
> + <entity index="70" value="871"/>
> + <entity index="71" value="873"/>
> + <entity index="72" value="874"/>
> + <entity index="73" value="875"/>
> + <entity index="74" value="876"/>
> + <entity index="75" value="877"/>
> + <entity index="76" value="878"/>
> + <entity index="77" value="879"/>
> + <entity index="78" value="880"/>
> + <entity index="79" value="881"/>
> + <entity index="80" value="882"/>
> + <entity index="81" value="883"/>
> + <entity index="82" value="884"/>
> + <entity index="83" value="885"/>
> + <entity index="84" value="886"/>
> + <entity index="85" value="887"/>
> + <entity index="86" value="888"/>
> + <entity index="87" value="889"/>
> + <entity index="88" value="890"/>
> + <entity index="89" value="891"/>
> + <entity index="90" value="922"/>
> + <entity index="91" value="921"/>
> + <entity index="92" value="923"/>
> + <entity index="93" value="963"/>
> + <entity index="94" value="962"/>
> + <entity index="95" value="964"/>
> + <entity index="96" value="1004"/>
> + <entity index="97" value="1003"/>
> + <entity index="98" value="1005"/>
> + <entity index="99" value="1045"/>
> + <entity index="100" value="1044"/>
> + <entity index="101" value="1046"/>
> + <entity index="102" value="1086"/>
> + <entity index="103" value="1085"/>
> + <entity index="104" value="1087"/>
> + <entity index="105" value="1127"/>
> + <entity index="106" value="1126"/>
> + <entity index="107" value="1128"/>
> + <entity index="108" value="1168"/>
> + <entity index="109" value="1167"/>
> + <entity index="110" value="1169"/>
> + <entity index="111" value="1209"/>
> + <entity index="112" value="1208"/>
> + <entity index="113" value="1210"/>
> + <entity index="114" value="1250"/>
> + <entity index="115" value="1249"/>
> + <entity index="116" value="1251"/>
> + </meshfunction>
> + </data_entry>
> + </data>
> + </mesh>
> +</dolfin>
>
> === added file 'nitsche_method/data/meshes/rotor_cavity.geo'
> --- nitsche_method/data/meshes/rotor_cavity.geo 1970-01-01 00:00:00 +0000
> +++ nitsche_method/data/meshes/rotor_cavity.geo 2010-03-16 14:39:09 +0000
> @@ -0,0 +1,82 @@
> +//This file defines start geometry for the driven rotor in the cavity problem.
> +
> +lb = 0.5;
> +lc = 0.3;
> +
> +// This variable can then be used in the definition of Gmsh's simplest
> +// `elementary entity', a `Point'. A Point is defined by a list of
> +// four numbers: three coordinates (X, Y and Z), and a characteristic
> +// length (lc) that sets the target element size at the point:
> +
> +//Cavity geometry
> +H = 2.0;
> +
> +Point(1) = {-H, -H, 0, lb};
> +Point(2) = {H, -H, 0, lb};
> +Point(3) = {H, H, 0, lb};
> +Point(4) = {-H, H, 0, lb};
> +
> +Line(1) = {1, 2};
> +Line(2) = {2, 3};
> +Line(3) = {3, 4};
> +Line(4) = {4, 1};
> +
> +Line Loop(5) = {1,2,3,4};
> +
> +//Physical Point(1) = {1,2,3,4};
> +//Physical Line(1) = {5};
> +
> +
> +//Rotor geometry
> +h = 0.1;
> +l = 1.0;
> +
> +//Point(5) = {-h,-h,0,lc};
> +//Point(6) = {h,-h,0,lc};
> +//Point(7) = {h,h,0,lc};
> +//Point(8) = {-h,h,0,lc};
> +//
> +//Line(6) = {5,6};
> +//Line(7) = {6,7};
> +//Line(8) = {7,8};
> +//Line(9) = {8,5};
> +//
> +//Line Loop(10) = {6,7,8,9};
> +//Plane Surface(101) = {10};
> +//Plane Surface(100) = {10};
> +//Plane Surface(200) = {5,10};
> +
> +Point(5) = {h,h,0,lc};
> +Point(6) = {h,l,0,lc};
> +Point(7) = {-h,l,0,lc};
> +Point(8) = {-h,h,0,lc};
> +Point(9) = {-l,h,0,lc};
> +Point(10) = {-l,-h,0,lc};
> +Point(11) = {-h,-h,0,lc};
> +Point(12) = {-h,-l,0,lc};
> +Point(13) = {h,-l,0,lc};
> +Point(14) = {h,-h,0,lc};
> +Point(15) = {l,-h,0,lc};
> +Point(16) = {l,h,0,lc};
> +
> +Line(6) = {5,6};
> +Line(7) = {6,7};
> +Line(8) = {7,8};
> +Line(9) = {8,9};
> +Line(10) = {9,10};
> +Line(11) = {10,11};
> +Line(12) = {11,12};
> +Line(13) = {12,13};
> +Line(14) = {13,14};
> +Line(15) = {14,15};
> +Line(16) = {14,15};
> +Line(17) = {15,16};
> +Line(18) = {16,5};
> +
> +Line Loop(19) = {6,7,8,9,10,11,12,13,14,15,16,17,18};
> +
> +//Physical Point(2) = {6,7,8,9,10,11,12,13,14,15,16,17,18};
> +//Physical Line(2) = {19};
> +
> +Plane Surface(100) = {19};
> +Plane Surface(101) = {5,19};
>
> === added directory 'nitsche_method/entity_intersection'
> === added file 'nitsche_method/entity_intersection/SConstruct'
> --- nitsche_method/entity_intersection/SConstruct 1970-01-01 00:00:00 +0000
> +++ nitsche_method/entity_intersection/SConstruct 2010-02-10 03:14:09 +0000
> @@ -0,0 +1,16 @@
> +# This is an example of a simple SConstruct file for building
> +# programs against DOLFIN. To build this demo, just type 'scons'.
> +
> +import os, commands
> +
> +# Get compiler from pkg-config
> +compiler = commands.getoutput('pkg-config --variable=compiler dolfin')
> +
> +# Create a SCons Environment based on the main os environment
> +env = Environment(ENV=os.environ, CXX=compiler)
> +
> +# Get compiler flags from pkg-config
> +env.ParseConfig('pkg-config --cflags --libs dolfin')
> +
> +#env.Program(['tetrahedron_triangle_intersection_test.cpp'])
> +env.Program(['primitive_intersector_test.cpp'])
>
> === added file 'nitsche_method/entity_intersection/primitive_intersector_test.cpp'
> --- nitsche_method/entity_intersection/primitive_intersector_test.cpp 1970-01-01 00:00:00 +0000
> +++ nitsche_method/entity_intersection/primitive_intersector_test.cpp 2010-03-05 23:36:16 +0000
> @@ -0,0 +1,40 @@
> +// =====================================================================================
> +//
> +// Copyright (C) 2010-02-10 André Massing
> +// Licensed under the GNU LGPL Version 2.1.
> +//
> +// Modified by André Massing, 2010
> +//
> +// First added: 2010-02-10
> +// Last changed: 2010-02-22
> +//
> +//Author: André Massing (am), massing@xxxxxxxxx
> +//Company: Simula Research Laboratory, Fornebu, Norway
> +//
> +// =====================================================================================
> +
> +#include <dolfin/mesh/dolfin_mesh.h>
> +#include <iostream>
> +
> +using namespace dolfin;
> +
> +int main ()
> +{
> + UnitCube cube(3,3,2);
> + cout <<"Total number of cells in Cube:" << cube.num_cells() <<endl;
> +
> + UnitSphere sphere(3);
> + cout <<"Total number of cells in Sphere:" << sphere.num_cells() <<endl;
> +
> + for (CellIterator cube_cell(cube); !cube_cell.end(); ++cube_cell)
> + {
> + for (CellIterator sphere_cell(cube); !sphere_cell.end(); ++sphere_cell)
> + {
> + if (PrimitiveIntersector::do_intersect(*cube_cell, *sphere_cell))
> + std::cout <<"Cell number " << cube_cell->index() << " intersects with "
> + <<"Cell number " << sphere_cell->index() << std::endl;
> + }
> + }
> +
> + return 0;
> +}
>
> === added directory 'nitsche_method/integration_on_complex_domains'
> === added directory 'nitsche_method/integration_on_complex_domains/gauss_projection_green'
> === added file 'nitsche_method/integration_on_complex_domains/gauss_projection_green/BaryCenterQuadrature.cpp'
> --- nitsche_method/integration_on_complex_domains/gauss_projection_green/BaryCenterQuadrature.cpp 1970-01-01 00:00:00 +0000
> +++ nitsche_method/integration_on_complex_domains/gauss_projection_green/BaryCenterQuadrature.cpp 2010-04-01 10:43:19 +0000
> @@ -0,0 +1,270 @@
> +// =====================================================================================
> +//
> +// Copyright (C) 2010-03-17 André Massing
> +// Licensed under the GNU LGPL Version 2.1.
> +//
> +// Modified by André Massing, 2010
> +//
> +// First added: 2010-03-17
> +// Last changed: 2010-04-01
> +//
> +//Author: André Massing (am), massing@xxxxxxxxx
> +//Company: Simula Research Laboratory, Fornebu, Norway
> +//Remark: This is an adapted version of Brian Mirtichs original C code for
> +//computin polyhedral mass properties, which can be found at <web-address>
> +// =====================================================================================
> +
> +#include "BaryCenterQuadrature.h"
> +#include <math.h>
> +
> +#define X 0
> +#define Y 1
> +#define Z 2
> +#define SQR(x) ((x)*(x))
> +#define CUBE(x) ((x)*(x)*(x))
> +
> +
> +using std::sqrt;
> +using namespace dolfin;
> +
> +
> +void BaryCenterQuadrature::compute_quadrature_rule(const Nef_polyhedron_3 & polyhedron)
> +{
> + //coordinates of the bary center
> + double bary_center [3] = {0,0,0};
> +
> + _weight = 0;
> +
> + //determine if polyhedron has a non vanishing 3d measure.
> + //This requires to have volumes marked with one.
> + bool has_volume = false;
> + Nef_polyhedron_3::Volume_const_iterator vi;
> + CGAL_forall_volumes(vi,polyhedron)
> + {
> + std::cout <<"vi->mark: " << vi->mark() <<std::endl;
> + if(vi->mark())
> + {
> + has_volume = true;
> + break;
> + }
> + }
> +
> + //variables needing during the computation
> + //components of surface normal
> + double nx, ny, nz;
> + double normal [3] = {0,0,0};
> +
> + //start and end point of the edge
> + double source[3] = {0,0,0};
> + double target[3] = {0,0,0};
> +
> + //permutation of the coordinates X,Y,Z
> + int A,B,C;
> +
> + std::cout <<" Number of volumes :" << polyhedron.number_of_volumes() << std::endl;
> + std::cout << "Outer iteration along faces, next f \n";
> + std::cout <<"++++++++++++++++++++++++++++++++++++++++" << std::endl;
> +
> + std::cout << "Polyhedron is 2d Manifold: " << const_cast<Nef_polyhedron_3 &>(polyhedron).is_simple() << std::endl;
> +
> + for(Nef_polyhedron_3::Halffacet_const_iterator f = polyhedron.halffacets_begin (),
> + end = polyhedron.halffacets_end();
> + f != end; ++f)
> + {
> + if (f->incident_volume()->mark() || (!has_volume && f->is_twin()))
> + {
> + //compute best projection direction.
> + //get normal of the surface and the offset
> + //FIXME: check right orientation
> + Plane_3 face_plane = f->plane();
> + Vector_3 test_normal = face_plane.orthogonal_vector();
> + nx = CGAL::to_double(face_plane.a());
> + ny = CGAL::to_double(face_plane.b());
> + nz = CGAL::to_double(face_plane.c());
> +
> + double normal_length = sqrt(nx*nx + ny*ny + nz*nz);
> + normal[0] = nx/normal_length;
> + normal[1] = ny/normal_length;
> + normal[2] = nz/normal_length;
> +
> + nx = fabs(nx);
> + ny = fabs(ny);
> + nz = fabs(nz);
> +
> + //Choose permutation to get best conditioned projection direction.
> + if (nx > ny && nx > nz)
> + C = X;
> + else
> + C = (ny > nz) ? Y : Z;
> + A = (C + 1) % 3;
> + B = (A + 1) % 3;
> +
> + //begin compute the face integrals
> + //begin compute the projection integrals
> + double a0, a1, da;
> + double b0, b1, db;
> + double a0_2, b0_2;
> + double a0_3, b0_3;
> + double a0_4, b0_4;
> + double a1_2, b1_2;
> + double C1, Ca, Caa, Cb, Cbb, Cab;
> + double Kab;
> + double P1, Pa, Pb, Paa, Pbb, Pab;
> + P1 = Pa = Pb = Paa = Pab = Pbb = 0.0;
> + double Fa, Fb, Fc, Faa, Fbb, Fcc;
> +
> + for(Nef_polyhedron_3::Halffacet_cycle_const_iterator fc = f->facet_cycles_begin();
> + fc != f->facet_cycles_end(); ++fc)
> + {
> + //exclude a SHalfloops (whatever this is)
> + if ( fc.is_shalfedge() )
> + {
> + Nef_polyhedron_3::SHalfedge_const_handle h = fc;
> + Nef_polyhedron_3::SHalfedge_around_facet_const_circulator hc(h);
> + Nef_polyhedron_3::SHalfedge_around_facet_const_circulator he(hc);
> + CGAL_For_all(hc,he)
> + {
> + std::cout <<"----------------------------------------" << std::endl;
> +
> + Nef_polyhedron_3::SVertex_const_handle v = hc->source();
> +
> + const Point_3 & source_point(v->source()->point());
> + const Point_3 & target_point(v->target()->point());
> + std::cout <<"Source point " << source_point << std::endl;
> + std::cout <<"Target point " << target_point <<std::endl;
> +
> + source[0] = CGAL::to_double(source_point.x());
> + source[1] = CGAL::to_double(source_point.y());
> + source[2] = CGAL::to_double(source_point.z());
> +
> + target[0] = CGAL::to_double(target_point.x());
> + target[1] = CGAL::to_double(target_point.y());
> + target[2] = CGAL::to_double(target_point.z());
> +
> + a0 = source[A];
> + b0 = source[B];
> + a1 = target[A];
> + b1 = target[B];
> +
> + da = a1 - a0;
> + db = b1 - b0;
> +
> + a0_2 = a0 * a0;
> + a0_3 = a0_2 * a0;
> + a0_4 = a0_3 * a0;
> + a1_2 = a1 * a1;
> + b0_2 = b0 * b0;
> + b0_3 = b0_2 * b0;
> + b0_4 = b0_3 * b0;
> + b1_2 = b1 * b1;
> +
> + C1 = a1 + a0;
> + Ca = a1*C1 + a0_2;
> + Caa = a1*Ca + a0_3;
> + Cb = b1*(b1 + b0) + b0_2;
> + Cbb = b1*Cb + b0_3;
> + Cab = 3*a1_2 + 2*a1*a0 + a0_2;
> + Kab = a1_2 + 2*a1*a0 + 3*a0_2;
> +
> + P1 += db*C1;
> + Pa += db*Ca;
> + Paa += db*Caa;
> + Pb += da*Cb;
> + Pbb += da*Cbb;
> + Pab += db*(b1*Cab + b0*Kab);
> + }
> + std::cout <<"----------------------------------------" << std::endl;
> + }
> + }
> +
> + P1 /= 2.0;
> + Pa /= 6.0;
> + Paa /= 12.0;
> + Pb /= -6.0;
> + Pbb /= -12.0;
> + Pab /= 24.0;
> +
> + info("P1 : %e",P1);
> + info("Pa : %e",Pa);
> + info("Paa : %e",Paa);
> + info("Pb : %e",Pb);
> + info("Pbb : %e",Pbb);
> + info("Pab : %e",Pab);
> +
> + double k1, k2, k3;
> + k1 = 1 / normal[C]; k2 = k1 * k1; k3 = k2 * k1;
> + double w = CGAL::to_double(face_plane.d())/normal_length ;
> + info("Normal[0] : %e ", normal[0]);
> + info("Normal[1] : %e ", normal[1]);
> + info("Normal[2] : %e ", normal[2]);
> + info("w : %e ", w);
> +
> + Fa = k1 * Pa;
> + Fb = k1 * Pb;
> + Fc = -k2 * (normal[A]*Pa + normal[B]*Pb + w*P1);
> +
> + //debug
> + Faa = k1 * Paa;
> + Fbb = k1 * Pbb;
> + Fcc = k3 * (SQR(normal[A])*Paa + 2*normal[A]*normal[B]*Pab + SQR(normal[B])*Pbb
> + + w*(2*(normal[A]*Pa + normal[B]*Pb) + w*P1));
> +
> + info("Fa : %e",Fa);
> + info("Fb : %e",Fb);
> + info("Fc : %e",Fc);
> + info("Faa : %e",Faa);
> + info("Fbb : %e",Fbb);
> + info("Fcc : %e",Fcc);
> + //end computation face integrals
> +
> + if (has_volume)
> + {
> + double contribution = (normal[X] * ((A == X) ? Fa : ((B == X) ? Fb : Fc)));
> + info("Contribution: %e", contribution);
> + info("Normal[X]: %e", normal[X]);
> + info("A: %d", A);
> + info("B: %d", B);
> + info("C: %d", C);
> + _weight += normal[X] * ((A == X) ? Fa : ((B == X) ? Fb : Fc));
> + info("Volume: %e", _weight);
> +
> + bary_center[A] += normal[A] * Faa;
> + bary_center[B] += normal[B] * Fbb;
> + bary_center[C] += normal[C] * Fcc;
> + }
> + else
> + {
> + info("No volume: P1 is %e",P1);
> + info("No volume: normal[C] is %e", normal[C]);
> + //does not work yet
> + //@todo does not work yet.
> + _weight += -P1/normal[C];
> + info("No volume: weight is %e", _weight);
> + info("Fa : %e",Fa);
> + info("Fb : %e",Fb);
> + info("Fc : %e",Fc);
> + bary_center[A] += -Fa;
> + bary_center[B] += -Fb;
> + bary_center[C] += -Fc;
> + }
> + }
> + }
> +
> + //applies also if polyhedron has no no volume, since
> + //we no iterating two times (each halffacet).
> + if (has_volume)
> + {
> + bary_center[X] /= (2*_weight);
> + bary_center[Y] /= (2*_weight);
> + bary_center[Z] /= (2*_weight);
> + }
> + else
> + {
> + bary_center[X] /= _weight;
> + bary_center[Y] /= _weight;
> + bary_center[Z] /= _weight;
> + }
> +
> + //use these later directly.
> + _point = Point(3,bary_center);
> +}
>
> === added file 'nitsche_method/integration_on_complex_domains/gauss_projection_green/BaryCenterQuadrature.h'
> --- nitsche_method/integration_on_complex_domains/gauss_projection_green/BaryCenterQuadrature.h 1970-01-01 00:00:00 +0000
> +++ nitsche_method/integration_on_complex_domains/gauss_projection_green/BaryCenterQuadrature.h 2010-04-01 10:43:19 +0000
> @@ -0,0 +1,69 @@
> +// =====================================================================================
> +//
> +// Copyright (C) 2010-03-17 André Massing
> +// Licensed under the GNU LGPL Version 2.1.
> +//
> +// Modified by André Massing, 2010
> +//
> +// First added: 2010-03-17
> +// Last changed: 2010-03-30
> +//
> +//Author: André Massing (am), massing@xxxxxxxxx
> +//Company: Simula Research Laboratory, Fornebu, Norway
> +//
> +// =====================================================================================
> +
> +#ifndef __BARYCENTERQUADRATURE_H
> +#define __BARYCENTERQUADRATURE_H
> +
> +#include <CGAL/Exact_predicates_exact_constructions_kernel.h>
> +
> +#include <CGAL/Nef_polyhedron_3.h>
> +#include <CGAL/Polyhedron_3.h>
> +
> +#include <dolfin/mesh/Point.h>
> +
> +namespace dolfin {
> +
> + typedef CGAL::Exact_predicates_exact_constructions_kernel Kernel;
> + typedef CGAL::Polyhedron_3<Kernel> Polyhedron_3;
> + typedef CGAL::Nef_polyhedron_3<Kernel> Nef_polyhedron_3;
> + typedef Nef_polyhedron_3::Point_3 Point_3;
> + typedef Nef_polyhedron_3::Plane_3 Plane_3;
> + typedef Nef_polyhedron_3::Vector_3 Vector_3;
> +
> + ///This class computes the barycenter of an arbitrary polyhedron or
> + ///polygon and therefore allows for bary center quadrature on complex
> + ///polyhedrons. Note: barycenter quadrature is exact for polynom deg <= 1.
> + class BaryCenterQuadrature {
> +
> + public:
> +
> + BaryCenterQuadrature() : _weight(0), _point(0,0,0), _length(1) {}
> +
> + ///@todo split the ugly implementation up into some subroutines and removes
> + //all the macro etc stuff. Make also working polyhedrons which are
> + //degenerated.
> + void compute_quadrature_rule(const Nef_polyhedron_3 & polyhedron);
> +
> + //return pointer
> + ///Gives the points for the last computed Polyhedron.
> + Point * points() {return &_point;}
> + ///Gives the weights for the last computed Polyhedron.
> + double * weights() {return &_weight; }
> +
> + ///Number of quadrature points/weights.
> + uint length() {return _length;}
> +
> + private:
> +
> + double _weight;
> + Point _point;
> + uint _length;
> +
> +
> + };
> +
> +}
> +
> +#endif // ----- #ifndef __BARYCENTERQUADRATURE_H -----
>
> === added file 'nitsche_method/integration_on_complex_domains/gauss_projection_green/README'
> --- nitsche_method/integration_on_complex_domains/gauss_projection_green/README 1970-01-01 00:00:00 +0000
> +++ nitsche_method/integration_on_complex_domains/gauss_projection_green/README 2010-03-18 18:15:26 +0000
> @@ -0,0 +1,146 @@
> +1. OVERVIEW
> +
> + This code accompanies the paper:
> +
> + Brian Mirtich, "Fast and Accurate Computation of
> + Polyhedral Mass Properties," journal of graphics
> + tools, volume 1, number 2, 1996.
> +
> + It computes the ten volume integrals needed for
> + determining the center of mass, moments of
> + inertia, and products of inertia for a uniform
> + density polyhedron. From this information, a
> + body frame can be computed.
> +
> + To compile the program, use an ANSI compiler, and
> + type something like
> +
> + % cc volInt.c -O2 -lm -o volInt
> +
> +
> + Revision history
> +
> + 26 Jan 1996 Program creation.
> +
> + 3 Aug 1996 Corrected bug arising when polyhedron density
> + is not 1.0. Changes confined to function main().
> + Thanks to Zoran Popovic for catching this one.
> +
> +
> +
> +2. POLYHEDRON GEOMETRY FILES
> +
> + The program reads a data file specified on the
> + command line. This data file describes the
> + geometry of a polyhedron, and has the following
> + format:
> +
> + N
> +
> + x_0 y_0 z_0
> + x_1 y_1 z_1
> + .
> + .
> + .
> + x_{N-1} y_{N-1} z_{N-1}
> +
> + M
> +
> + k1 v_{1,1} v_{1,2} ... v_{1,k1}
> + k2 v_{2,1} v_{2,2} ... v_{2,k2}
> + .
> + .
> + .
> + kM v_{M,1} v_{M,2} ... v_{M,kM}
> +
> +
> + where:
> + N number of vertices on polyhedron
> + x_i y_i z_i x, y, and z coordinates of ith vertex
> + M number of faces on polyhedron
> + ki number of vertices on ith face
> + v_{i,j} jth vertex on ith face
> +
> + In English:
> +
> + First the number of vertices are specified. Next
> + the vertices are defined by listing the 3
> + coordinates of each one. Next the number of faces
> + are specified. Finally, the faces themselves are
> + specified. A face is specified by first giving
> + the number of vertices around the polygonal face,
> + followed by the (integer) indices of these
> + vertices. When specifying indices, note that
> + they must be given in counter-clockwise order
> + (when looking at the face from outside the
> + polyhedron), and the vertices are indexed from 0
> + to N-1 for a polyhedron with N faces.
> +
> + White space can be inserted (or not) as desired.
> + Three example polyhedron geometry files are included:
> +
> + cube A cube, 20 units on a side, centered at
> + the origin and aligned with the coordinate axes.
> +
> + tetra A tetrahedron formed by taking the convex
> + hull of the origin, and the points (5,0,0),
> + (0,4,0), and (0,0,3).
> +
> + icosa An icosahedron with vertices lying on the unit
> + sphere, centered at the origin.
> +
> +
> +
> +3. RUNNING THE PROGRAM
> +
> + Simply type,
> +
> + % volInt <polyhedron geometry filename>
> +
> + The program will read in the geometry of the
> + polyhedron, and the print out the ten volume
> + integrals.
> +
> + The program also computes some of the mass
> + properties which may be inferred from the volume
> + integrals. A density of 1.0 is assumed, although
> + this may be changed in function main(). The
> + center of mass is computed as well as the inertia
> + tensor relative to a frame with origin at the
> + center of mass. Note, however, that this matrix
> + may not be diagonal. If not, a diagonalization
> + routine must be performed, to determine the
> + orientation of the true body frame relative to
> + the original model frame. The Jacobi method
> + works quite well (see Numerical Recipes in C by
> + Press, et. al.).
> +
> +
> +
> +4. DISCLAIMERS
> +
> + 1. The volume integration code has been written
> + to match the development and algorithms presented
> + in the paper, and not with maximum optimization
> + in mind. While inherently very efficient, a few
> + more cycles can be squeezed out of the algorithm.
> + This is left as an exercise. :)
> +
> + 2. Don't like global variables? The three
> + procedures which evaluate the volume integrals
> + can be combined into a single procedure with two
> + nested loops. In addition to providing some
> + speedup, all of the global variables can then be
> + made local.
> +
> + 3. The polyhedron data structure used by the
> + program is admittedly lame; much better schemes
> + are possible. The idea here is just to give the
> + basic integral evaluation code, which will have
> + to be adjusted for other polyhedron data
> + structures.
> +
> + 4. There is no error checking for the input
> + files. Be careful. Note the hard limits
> + #defined for the number of vertices, number of
> + faces, and number of vertices per faces.
>
> === added file 'nitsche_method/integration_on_complex_domains/gauss_projection_green/SConstruct'
> --- nitsche_method/integration_on_complex_domains/gauss_projection_green/SConstruct 1970-01-01 00:00:00 +0000
> +++ nitsche_method/integration_on_complex_domains/gauss_projection_green/SConstruct 2010-04-02 18:46:57 +0000
> @@ -0,0 +1,22 @@
> +import os, commands
> +
> +# Import("cppunitCfg")
> +
> +# Get compiler from pkg-config
> +compiler = commands.getoutput('pkg-config --variable=compiler dolfin')
> +
> +# Create a SCons Environment based on the main os environment
> +env = Environment(ENV=os.environ, CXX=compiler)
> +
> +env.Append(CXXFLAGS = '-ggdb')
> +
> +# Get compiler flags from pkg-config
> +env.ParseConfig('pkg-config --cflags --libs dolfin')
> +
> +#Program
> +env.Program('quadrature_test',['BaryCenterQuadrature.cpp','main.cpp'])
> +env.Program('iteration_tetrahedron.cpp')
> +env.Program('iteration_tetrahedron_2.cpp')
> +
> +env.Append(LIBS='-lcppunit')
> +env.Program('unittest',['BaryCenterQuadrature.cpp','test.cpp'])
>
> === added file 'nitsche_method/integration_on_complex_domains/gauss_projection_green/alut_snippet.cpp'
> --- nitsche_method/integration_on_complex_domains/gauss_projection_green/alut_snippet.cpp 1970-01-01 00:00:00 +0000
> +++ nitsche_method/integration_on_complex_domains/gauss_projection_green/alut_snippet.cpp 2010-03-18 18:15:26 +0000
> @@ -0,0 +1,96 @@
> +// =====================================================================================
> +//
> +// Copyright (C) 2010-02-17 André Massing
> +// Licensed under the GNU LGPL Version 2.1.
> +//
> +// Modified by André Massing, 2010
> +//
> +// First added: 2010-02-17
> +// Last changed: 2010-02-17
> +//
> +//Author: André Massing (am), massing@xxxxxxxxx
> +//Company: Simula Research Laboratory, Fornebu, Norway
> +//
> +// =====================================================================================
> +
> +
> +struct
> +{
> + Point_3 p1, p2, p3;
> +}triangle;
> +
> +class Shell_explorer
> +{
> + bool first;
> + const Nef_polyhedron& N;
> + vector<triangle> m_triangle_list;
> +
> +public:
> + Shell_explorer(const Nef_polyhedron& N_)
> + : first(true), N(N_) {}
> +
> + void visit(Vertex_const_handle v){}
> + void visit(Halfedge_const_handle ){}
> + void visit(Halffacet_const_handle facet)
> + {
> + for (Halffacet_cycle_const_iterator it =
> +facet->facet_cycles_begin(); it != facet->facet_cycles_end(); it++)
> + {
> + SHalfedge_const_handle se = SHalfedge_const_handle(it);
> + CGAL_assertion(se!=0);
> + SHalfedge_around_facet_const_circulator hc_start(se);
> + SHalfedge_around_facet_const_circulator hc_end(hc_start);
> + vector <Point_3> vertices;
> + int count = 0;
> + CGAL_For_all(hc_start,hc_end)
> + {
> + SVertex_const_handle svert = hc_start->source();
> + Point_3 vpoint = svert->center_vertex()->point();
> + vertices.push_back(vpoint);
> + }
> + for (int i=0; i<vertices.size()-2; i++)
> + {
> + triangle t;
> + t.p1 = vertices[0];
> + t.p3 = vertices[i+1];
> + t.p2 = vertices[i+2];
> + m_triangle_list.push_back(t);
> + }
> +
> + }
> + }
> + void visit(SHalfedge_const_handle ) {}
> + void visit(SHalfloop_const_handle ) {}
> + void visit(SFace_const_handle sf) {}
> +
> + vector<triangle>getTriangleList()
> + {
> + vector<triangle> templist;
> + for (int i=0; i<m_triangle_list.size(); i++)
> + templist.push_back(m_triangle_list[i]);
> + m_triangle_list.clear();
> + return templist;
> + }
> + void reset_minimal_vertex() { first = true; }
> +};
> +
> +Usage:
> +void main()
> +{
> + Nef_polyhedron nPoly;
> + Volume_const_iterator c;
> + Shell_explorer SE(nPoly);
> +
> + CGAL_forall_volumes(c,nPoly)
> + {
> + Shell_entry_const_iterator it;
> + int shellcount = 0;
> + vector <triangle> temp;
> + CGAL_forall_shells_of(it,c)
> + {
> + SE.reset_minimal_vertex();
> + entrypaths.visit_shell_objects(SFace_const_handle(it),SE);
> + vector<triangle> tris = SE.getTriangleList();
> + }
> + }
> +}
>
> === added file 'nitsche_method/integration_on_complex_domains/gauss_projection_green/cube'
> --- nitsche_method/integration_on_complex_domains/gauss_projection_green/cube 1970-01-01 00:00:00 +0000
> +++ nitsche_method/integration_on_complex_domains/gauss_projection_green/cube 2010-03-18 18:15:26 +0000
> @@ -0,0 +1,21 @@
> +8
> +
> +-10 -10 -10
> ++10 -10 -10
> ++10 +10 -10
> +-10 +10 -10
> +-10 -10 +10
> ++10 -10 +10
> ++10 +10 +10
> +-10 +10 +10
> +
> +6
> +
> +4 0 3 2 1
> +4 4 5 6 7
> +4 0 1 5 4
> +4 6 2 3 7
> +4 1 2 6 5
> +4 0 4 7 3
> +
> +
>
> === added file 'nitsche_method/integration_on_complex_domains/gauss_projection_green/icosa'
> --- nitsche_method/integration_on_complex_domains/gauss_projection_green/icosa 1970-01-01 00:00:00 +0000
> +++ nitsche_method/integration_on_complex_domains/gauss_projection_green/icosa 2010-03-18 18:15:26 +0000
> @@ -0,0 +1,38 @@
> +12
> +
> ++0.00000000 +0.00000000 +1.00000000
> ++0.00000000 +0.00000000 -1.00000000
> ++0.89442719 +0.00000000 +0.44721360
> ++0.27639320 +0.85065081 +0.44721360
> +-0.72360680 +0.52573111 +0.44721360
> +-0.72360680 -0.52573111 +0.44721360
> ++0.27639320 -0.85065081 +0.44721360
> ++0.72360680 +0.52573111 -0.44721360
> +-0.27639320 +0.85065081 -0.44721360
> +-0.89442719 +0.00000000 -0.44721360
> +-0.27639320 -0.85065081 -0.44721360
> ++0.72360680 -0.52573111 -0.44721360
> +
> +20
> +
> +3 6 11 2
> +3 3 2 7
> +3 7 2 11
> +3 0 2 3
> +3 0 6 2
> +3 10 1 11
> +3 1 7 11
> +3 10 11 6
> +3 8 7 1
> +3 8 3 7
> +3 5 10 6
> +3 5 6 0
> +3 4 3 8
> +3 4 0 3
> +3 9 8 1
> +3 9 1 10
> +3 4 5 0
> +3 9 10 5
> +3 9 5 4
> +3 9 4 8
> +
>
> === added file 'nitsche_method/integration_on_complex_domains/gauss_projection_green/iteration_tetrahedron.cpp'
> --- nitsche_method/integration_on_complex_domains/gauss_projection_green/iteration_tetrahedron.cpp 1970-01-01 00:00:00 +0000
> +++ nitsche_method/integration_on_complex_domains/gauss_projection_green/iteration_tetrahedron.cpp 2010-03-25 23:33:08 +0000
> @@ -0,0 +1,58 @@
> +#include <CGAL/Exact_predicates_exact_constructions_kernel.h>
> +#include <CGAL/Nef_polyhedron_3.h>
> +#include <CGAL/Polyhedron_3.h>
> +
> +typedef CGAL::Exact_predicates_exact_constructions_kernel Kernel;
> +typedef CGAL::Polyhedron_3<Kernel> Polyhedron_3;
> +typedef CGAL::Nef_polyhedron_3<Kernel> Nef_polyhedron_3;
> +typedef Nef_polyhedron_3::Point_3 Point_3;
> +
> +int main()
> +{
> +
> + Point_3 a(1.0, 0.0, 0.0);
> + Point_3 b(0.0, 1.0, 0.0);
> + Point_3 c(0.0, 0.0, 1.0);
> + Point_3 d(0.0, 0.0, 0.0);
> +
> + //Build a tetrahedron
> + Polyhedron_3 P;
> + P.make_tetrahedron(a,b,c,d);
> + Nef_polyhedron_3 NP(P);
> +
> + std::cout << "Iteration along face ...\n";
> + int face_counter = 0;
> + for(Nef_polyhedron_3::Halffacet_const_iterator f = NP.halffacets_begin (),
> + end = NP.halffacets_end();
> + f != end;
> + ++f)
> + {
> + if(f->is_twin())
> + continue;
> + std::cout <<"########################################" << std::endl;
> + std::cout << "Face number " << face_counter++ << std::endl;
> + for(Nef_polyhedron_3::Halffacet_cycle_const_iterator fc = f->facet_cycles_begin(),
> + cycles_end = f->facet_cycles_end();
> + fc != cycles_end;
> + ++fc)
> + {
> + if ( fc.is_shalfedge() )
> + {
> + std::cout <<"Edge consist of points: " << std::endl;
> + Nef_polyhedron_3::SHalfedge_const_handle se = fc;
> + Nef_polyhedron_3::SHalfedge_around_facet_const_circulator hc(se);
> + Nef_polyhedron_3::SHalfedge_around_facet_const_circulator he(hc);
> + CGAL_For_all(hc,he)
> + {
> + Nef_polyhedron_3::SVertex_const_handle v = hc->source();
> +// Nef_polyhedron_3::Halfedge_const_handle v = hc->source();
> + const Point_3 & source_point(v->source()->point());
> + const Point_3 & target_point(v->target()->point());
> + std::cout <<"Source point " << source_point << std::endl;
> + std::cout <<"Target point " << target_point <<std::endl;
> + }
> + }
> + }
> + }
> + return 0;
> +}
>
> === added file 'nitsche_method/integration_on_complex_domains/gauss_projection_green/iteration_tetrahedron_2.cpp'
> --- nitsche_method/integration_on_complex_domains/gauss_projection_green/iteration_tetrahedron_2.cpp 1970-01-01 00:00:00 +0000
> +++ nitsche_method/integration_on_complex_domains/gauss_projection_green/iteration_tetrahedron_2.cpp 2010-04-01 10:43:19 +0000
> @@ -0,0 +1,140 @@
> +#include<vector>
> +
> +#include <CGAL/Exact_predicates_exact_constructions_kernel.h>
> +#include <CGAL/Nef_polyhedron_3.h>
> +#include <CGAL/Polyhedron_3.h>
> +#include <CGAL/Gmpz.h>
> +#include <CGAL/Homogeneous.h>
> +
> +typedef CGAL::Exact_predicates_exact_constructions_kernel Kernel;
> +typedef CGAL::Polyhedron_3<Kernel> Polyhedron_3;
> +typedef CGAL::Nef_polyhedron_3<Kernel> Nef_polyhedron_3;
> + typedef Nef_polyhedron_3::Vector_3 Vector_3;
> +typedef Nef_polyhedron_3::Point_3 Point_3;
> +typedef Nef_polyhedron_3::Vertex_const_handle Vertex_const_handle;
> +typedef Nef_polyhedron_3::Halfedge_const_handle Halfedge_const_handle;
> +typedef Nef_polyhedron_3::Halffacet_const_handle Halffacet_const_handle;
> +typedef Nef_polyhedron_3::Halffacet_cycle_const_iterator Halffacet_cycle_const_iterator;
> +typedef Nef_polyhedron_3::SHalfedge_const_handle SHalfedge_const_handle;
> +typedef Nef_polyhedron_3::SHalfedge_around_facet_const_circulator SHalfedge_around_facet_const_circulator;
> +typedef Nef_polyhedron_3::SHalfloop_const_handle SHalfloop_const_handle;
> +typedef Nef_polyhedron_3::SFace_const_handle SFace_const_handle;
> +typedef Nef_polyhedron_3::Volume_const_iterator Volume_const_iterator;
> +typedef Nef_polyhedron_3::Shell_entry_const_iterator Shell_entry_const_iterator;
> +
> +typedef std::vector<Nef_polyhedron_3> PolyhedronList;
> +typedef PolyhedronList::const_iterator PolyhedronListIterator;
> +
> +class Shell_explorer
> +{
> + const Nef_polyhedron_3 & N;
> +
> +public:
> + Shell_explorer(const Nef_polyhedron_3 & N_)
> + : N(N_) {}
> +
> + void visit(Vertex_const_handle v) {}
> + void visit(Halfedge_const_handle ){}
> + void visit(Halffacet_const_handle facet)
> + {
> + std::cout << "Visting faces ..." << std::endl;
> + std::cout <<"++++++++++++++++++++++++++++++++++++++++" << std::endl;
> + Vector_3 test_normal = facet->plane().orthogonal_vector();
> + std::cout <<"Normal vector of the plane is " << test_normal <<std::endl;
> +
> + for (Halffacet_cycle_const_iterator it =
> + facet->facet_cycles_begin(); it != facet->facet_cycles_end(); it++)
> + {
> + std::cout << "Visting Edges ..." << std::endl;
> + std::cout <<"----------------------------------------" << std::endl;
> + Nef_polyhedron_3::SHalfedge_const_handle h = it;
> + Nef_polyhedron_3::SHalfedge_around_facet_const_circulator hc(h), he(hc);
> + CGAL_For_all(hc,he)
> + {
> + Nef_polyhedron_3::SVertex_const_handle v = hc->source();
> + const Nef_polyhedron_3::Point_3& source = v->source()->point();
> + const Nef_polyhedron_3::Point_3& target = v->target()->point();
> + std::cout <<"Source :" <<source <<std::endl;
> + std::cout <<"Target :" <<target <<std::endl;
> + }
> + }
> + }
> + void visit(SHalfedge_const_handle ) {}
> + void visit(SHalfloop_const_handle ) {}
> + void visit(SFace_const_handle sf) {}
> +
> +};
> +
> +int main()
> +{
> + Point_3 a1(1.0, 0.0, 0.0);
> + Point_3 b1(0.0, 1.0, 0.0);
> + Point_3 c1(0.0, 0.0, 1.0);
> + Point_3 d1(0.0, 0.0, 0.0);
> +
> + Point_3 a2(-1.0, 0.0, 0.0);
> + Point_3 b2(0.0, 1.0, 0.0);
> + Point_3 c2(0.0, 0.0, 1.0);
> + Point_3 d2(0.0, 0.0, 0.0);
> +
> + Point_3 a3(1.0, 0.0, 0.0);
> + Point_3 b3(0.0, -1.0, 0.0);
> + Point_3 c3(0.0, 0.0, 1.0);
> + Point_3 d3(0.0, 0.0, 0.0);
> +
> + Point_3 a4(-1.0, 0.0, 0.0);
> + Point_3 b4(0.0, -1.0, 0.0);
> + Point_3 c4(0.0, 0.0, 1.0);
> + Point_3 d4(0.0, 0.0, 0.0);
> +
> + Point_3 A(2.0, 1.0, 1.0);
> + Point_3 B(1.0, 2.0, 1.0);
> + Point_3 C(1.0, 1.0, 2.0);
> + Point_3 D(1.0, 1.0, 1.0);
> +
> + //Build a point tetrahedron
> + Polyhedron_3 P1;
> + P1.make_triangle(a1,b1,c1);
> + Nef_polyhedron_3 N1(P1);
> +
> + Polyhedron_3 P2;
> + P2.make_tetrahedron(a2,b2,c2,d2);
> + Nef_polyhedron_3 N2(P2);
> +
> + Polyhedron_3 P3;
> + P3.make_tetrahedron(a3,b3,c3,d3);
> + Nef_polyhedron_3 N3(P3);
> +
> + Polyhedron_3 P4;
> + P4.make_tetrahedron(a4,b4,c4,d4);
> + Nef_polyhedron_3 N4(P4);
> +
> + PolyhedronList polyhedrons;
> +
> + polyhedrons.push_back(N1);
> + polyhedrons.push_back(N2);
> +
> + //Build a point tetrahedron
> + for (PolyhedronListIterator pi = polyhedrons.begin(), end = polyhedrons.end();
> + pi != end; ++pi)
> + {
> + Volume_const_iterator vi;
> + Shell_explorer SE(*pi);
> +
> + std::cout <<"Next polyhedron:\n"
> + << "##############################" <<std::endl;
> + int ic = 0;
> + CGAL_forall_volumes(vi,*pi)
> + {
> + std::cout << "Volume " << ic++ << " has mark " << vi->mark()<< std::endl;
> +
> + Shell_entry_const_iterator it;
> + CGAL_forall_shells_of(it,vi)
> + {
> + pi->visit_shell_objects(SFace_const_handle(it),SE);
> + }
> + std::cout << "\n\n";
> + }
> + }
> + return 0;
> +}
>
> === added file 'nitsche_method/integration_on_complex_domains/gauss_projection_green/main.cpp'
> --- nitsche_method/integration_on_complex_domains/gauss_projection_green/main.cpp 1970-01-01 00:00:00 +0000
> +++ nitsche_method/integration_on_complex_domains/gauss_projection_green/main.cpp 2010-04-01 10:43:19 +0000
> @@ -0,0 +1,190 @@
> +// =====================================================================================
> +//
> +// Copyright (C) 2010-03-19 André Massing
> +// Licensed under the GNU LGPL Version 2.1.
> +//
> +// Modified by André Massing, 2010
> +//
> +// First added: 2010-03-19
> +// Last changed: 2010-03-31
> +//
> +//Author: André Massing (am), massing@xxxxxxxxx
> +//Company: Simula Research Laboratory, Fornebu, Norway
> +//
> +// =====================================================================================
> +
> +#include "BaryCenterQuadrature.h"
> +#include <dolfin.h>
> +
> +#include <iostream>
> +#include<vector>
> +
> +using namespace dolfin;
> +
> +typedef std::vector<Nef_polyhedron_3> PolyhedronList;
> +typedef PolyhedronList::const_iterator PolyhedronListIterator;
> +
> +
> +using namespace dolfin;
> +
> +int main()
> +{
> +
> + Point_3 e0(0.0, 0.0, 0.0);
> +
> + Point_3 e1(1.0, 0.0, 0.0);
> + Point_3 e2(0.0, 1.0, 0.0);
> + Point_3 e3(0.0, 0.0, 1.0);
> +
> + Point_3 _e1(-1.0, 0.0, 0.0);
> + Point_3 _e2(0.0, -1.0, 0.0);
> + Point_3 _e3(0.0, 0.0, -1.0);
> +
> + Point_3 A(2.0, 1.0, 1.0);
> + Point_3 B(1.0, 2.0, 1.0);
> + Point_3 C(1.0, 1.0, 2.0);
> + Point_3 D(1.0, 1.0, 1.0);
> +
> + //Create tetrahedrons with unit vectors
> + Polyhedron_3 P1;
> + P1.make_tetrahedron(e0,e1,e2,e3);
> + Nef_polyhedron_3 N1(P1);
> +
> + Polyhedron_3 P2;
> + P2.make_tetrahedron(e0,_e1,e2,e3);
> + Nef_polyhedron_3 N2(P2);
> +
> + Polyhedron_3 P3;
> + P3.make_tetrahedron(e0,e1,_e2,e3);
> + Nef_polyhedron_3 N3(P3);
> +
> + Polyhedron_3 P4;
> + P4.make_tetrahedron(e0,_e1,_e2,e3);
> + Nef_polyhedron_3 N4(P4);
> +
> + Polyhedron_3 P5;
> + P5.make_tetrahedron(e0,e1,e2,_e3);
> + Nef_polyhedron_3 N5(P5);
> +
> + Polyhedron_3 P6;
> + P6.make_tetrahedron(e0,_e1,e2,_e3);
> + Nef_polyhedron_3 N6(P6);
> +
> + Polyhedron_3 P7;
> + P7.make_tetrahedron(e0,e1,_e2,_e3);
> + Nef_polyhedron_3 N7(P7);
> +
> + Polyhedron_3 P8;
> + P8.make_tetrahedron(e0,_e1,_e2,_e3);
> + Nef_polyhedron_3 N8(P8);
> +
> + //Create triangles
> + Polyhedron_3 P9;
> + P9.make_triangle(e0,e1,e2);
> + Nef_polyhedron_3 N9(P9);
> +
> + Polyhedron_3 P10;
> + P10.make_triangle(e0,e1,e3);
> + Nef_polyhedron_3 N10(P10);
> +
> + Polyhedron_3 P11;
> + P11.make_triangle(e0,e2,e3);
> + Nef_polyhedron_3 N11(P11);
> +
> + Polyhedron_3 P12;
> + P12.make_triangle(e1,e2,e3);
> + Nef_polyhedron_3 N12(P12);
> +
> + Polyhedron_3 P13;
> + P13.make_triangle(e0,_e1,e2);
> + Nef_polyhedron_3 N13(P13);
> +
> + Polyhedron_3 P14;
> + P14.make_triangle(e0,_e1,e3);
> + Nef_polyhedron_3 N14(P14);
> +
> + Polyhedron_3 P15;
> + P15.make_triangle(e0,e2,e3);
> + Nef_polyhedron_3 N15(P15);
> +
> + Polyhedron_3 P16;
> + P16.make_triangle(_e1,e2,e3);
> + Nef_polyhedron_3 N16(P16);
> +
> + PolyhedronList polyhedrons;
> + PolyhedronList polygons;
> +
> + //single polyhedrons
> + polyhedrons.push_back(N1);
> + polyhedrons.push_back(N2);
> + polyhedrons.push_back(N3);
> + polyhedrons.push_back(N4);
> + polyhedrons.push_back(N5);
> + polyhedrons.push_back(N6);
> + polyhedrons.push_back(N7);
> + polyhedrons.push_back(N8);
> +
> + //combined tetrehedrons
> + polyhedrons.push_back(N1 + N2);
> + polyhedrons.push_back(N1 + N2 + N3);
> + polyhedrons.push_back(N1 + N2 + N3 + N4);
> +
> + polyhedrons.push_back(N5 + N6);
> + polyhedrons.push_back(N5 + N6 + N7);
> + polyhedrons.push_back(N5 + N6 + N7 + N9);
> +
> + polyhedrons.push_back(N1 + N2 + N3 + N4 + N5 + N6 + N7 + N4);
> +
> + //polygons
> + polygons.push_back(N9);
> + polygons.push_back(N10);
> + polygons.push_back(N11);
> + polygons.push_back(N12);
> + polygons.push_back(N13);
> + polygons.push_back(N14);
> + polygons.push_back(N15);
> + polygons.push_back(N16);
> +
> + //combined polygons
> + polygons.push_back(N9 + N10);
> + polygons.push_back(N9 + N10 + N11);
> + polygons.push_back(N9 + N10 + N11 + N12);
> +
> + polygons.push_back(N13 + N14);
> + polygons.push_back(N13 + N14 + N15);
> + polygons.push_back(N13 + N14 + N15 + N16);
> +
> + polygons.push_back(N9 + N10 + N11 + N12 + N13 + N14 + N15 + N12);
> +
> + BaryCenterQuadrature quadrature_rule;
> +
> + int counter = 0;
> + std::cout << "Quadrature rule for Polyhdrons\n" ;
> + std::cout <<"------------------------------\n";
> + for (PolyhedronListIterator pi = polyhedrons.begin(), end = polyhedrons.end();
> + pi != end; ++pi)
> + {
> + cout <<"Next polyhedron: " << ++counter << endl
> + << "##############################" <<endl;
> + quadrature_rule.compute_quadrature_rule(*pi);
> + cout <<"Volume: " << *(quadrature_rule.weights()) <<endl;
> + cout <<"Barycenter: "
> + << quadrature_rule.points()[0] << endl <<endl;
> + }
> +
> + counter = 0;
> + std::cout << "Quadrature rule for Polygons\n" ;
> + std::cout <<"------------------------------\n";
> + for (PolyhedronListIterator pi = polygons.begin(), end = polygons.end();
> + pi != end; ++pi)
> + {
> + cout <<"Next polygons: " << ++counter << endl
> + << "##############################" <<endl;
> + quadrature_rule.compute_quadrature_rule(*pi);
> + cout <<"Volume: " << *(quadrature_rule.weights()) <<endl;
> + cout <<"Barycenter: "
> + << quadrature_rule.points()[0] << endl <<endl;
> + }
> +
> + return 0;
> +}
>
> === added file 'nitsche_method/integration_on_complex_domains/gauss_projection_green/martin_baeker_snippet.cpp'
> --- nitsche_method/integration_on_complex_domains/gauss_projection_green/martin_baeker_snippet.cpp 1970-01-01 00:00:00 +0000
> +++ nitsche_method/integration_on_complex_domains/gauss_projection_green/martin_baeker_snippet.cpp 2010-03-18 18:15:26 +0000
> @@ -0,0 +1,93 @@
> +// =====================================================================================
> +//
> +// Copyright (C) 2010-02-17 André Massing
> +// Licensed under the GNU LGPL Version 2.1.
> +//
> +// Modified by André Massing, 2010
> +//
> +// First added: 2010-02-17
> +// Last changed: 2010-02-17
> +//
> +//Author: André Massing (am), massing@xxxxxxxxx
> +//Company: Simula Research Laboratory, Fornebu, Norway
> +//
> +// =====================================================================================
> +
> + for(Nef_polyhedron::Halffacet_const_iterator
> + f = Cube.halffacets_begin (),
> + end = Cube.halffacets_end();
> + f != end; ++f)
> + {
> + if(f->is_twin()) continue;
> + std::cout << "Outer iteration, next f \n";
> + for(Nef_polyhedron::Halffacet_cycle_const_iterator
> + fc = f->facet_cycles_begin(),
> + end = f->facet_cycles_end();
> + fc != end; ++fc)
> + {
> + std::cout << "Inner iteration, next fc \n";
> + if ( fc.is_shalfedge() )
> + {
> + std::cout << "Halfedge consists of points \n";
> +
> + Nef_polyhedron::SHalfedge_const_handle h = fc;
> + Nef_polyhedron::SHalfedge_around_facet_const_circulator hc(h), he(hc);
> + CGAL_For_all(hc,he)
> + { // all vertex coordinates in facet cycle
> + Nef_polyhedron::SVertex_const_handle v = hc->source();
> + const Nef_polyhedron::Point_3& point = v->source()->point();
> + std::cout << "p: " << CGAL::to_double(point.x()) << " " <<
> + CGAL::to_double(point.y()) << " " <<
> + CGAL::to_double(point.z()) << std::endl;
> +
> + }
> +
> + }
> + }
> +
> + }
> +
> + // Alternative to iterating over halfedges:
> +// Nef_polyhedron::Halfedge_const_iterator e = Cube.halfedges_begin();
> +// CGAL_forall_halfedges(e,Cube) {
> +
> +
> + for(Nef_polyhedron::Halfedge_const_iterator
> + e = Cube.halfedges_begin(),
> + end = Cube.halfedges_end();
> + e != end; ++e)
> + {
> + const Nef_polyhedron::Vertex_const_handle& s = e->source();
> + const Nef_polyhedron::Vertex_const_handle& t = e->twin()->source();
> + const Nef_polyhedron::Point_3& a = s->point();
> + const Nef_polyhedron::Point_3& b = t->point();
> + std::cout << "From a: " << CGAL::to_double(a.x()) << " " <<
> + CGAL::to_double(a.y()) << " " <<
> + CGAL::to_double(a.z()) << " to b: " <<
> + CGAL::to_double(b.x()) << " " <<
> + CGAL::to_double(b.y()) << " " <<
> + CGAL::to_double(b.z()) <<std::endl;
> + typedef Nef_polyhedron::Object_handle Object_handle;
> + typedef Nef_polyhedron::Vertex_const_handle Vertex_const_handle;
> + Vertex_const_handle v;
> + Object_handle o = Cube.locate(a);
> + if(CGAL::assign(v,o))
> + std::cout << " a is a vertex" << std::endl;
> + else
> + std::cout << " a is not a vertex" << std::endl;
> + }
> +
> +// Nef_polyhedron::Vertex_const_iterator v;
> +
> +// CGAL_forall_vertices(v, Cube.sncp())
> + for(Nef_polyhedron::Vertex_const_iterator
> + v = Cube.vertices_begin(),
> + end = Cube.vertices_end();
> + v != end; ++v)
> + {
> + const Nef_polyhedron::Point_3& a = v->point();
> + std::cout << "Vertex: " << CGAL::to_double(a.x()) << " " <<
> + CGAL::to_double(a.y()) << " " <<
> + CGAL::to_double(a.z()) << std::endl;
> +
> + }
>
> === added file 'nitsche_method/integration_on_complex_domains/gauss_projection_green/test.cpp'
> --- nitsche_method/integration_on_complex_domains/gauss_projection_green/test.cpp 1970-01-01 00:00:00 +0000
> +++ nitsche_method/integration_on_complex_domains/gauss_projection_green/test.cpp 2010-04-02 22:47:40 +0000
> @@ -0,0 +1,598 @@
> +// =====================================================================================
> +//
> +// Copyright (C) 2010-04-01 André Massing
> +// Licensed under the GNU LGPL Version 2.1.
> +//
> +// Modified by André Massing, 2010
> +//
> +// First added: 2010-04-01
> +// Last changed: 2010-04-03
> +//
> +//Author: André Massing (am), massing@xxxxxxxxx
> +//Company: Simula Research Laboratory, Fornebu, Norway
> +//
> +//Description: Unittest for BaryCenterQuadrature. =====================================================================================
> +
> +#include<vector>
> +#include <iostream>
> +
> +#include <dolfin.h>
> +#include <dolfin/common/unittest.h>
> +
> +#include <CGAL/Nef_polyhedron_3.h>
> +#include "BaryCenterQuadrature.h"
> +
> +using namespace dolfin;
> +
> +typedef Nef_polyhedron_3::Aff_transformation_3 Aff_transformation_3;
> +
> +typedef std::vector<Nef_polyhedron_3> PolyhedronList;
> +typedef PolyhedronList::const_iterator PolyhedronListIterator;
> +
> +typedef std::vector<int> IntList;
> +typedef std::vector<int>::const_iterator IntListIterator;
> +typedef std::vector<double> DoubleList;
> +typedef std::vector<double>::const_iterator DoubleListIterator;
> +typedef std::vector<Point> PointList;
> +typedef std::vector<Point>::const_iterator PointListIterator;
> +
> +class BaryCenter : public CppUnit::TestFixture
> +{
> + CPPUNIT_TEST_SUITE(BaryCenter);
> + CPPUNIT_TEST(testSimplePolyhedrons);
> + CPPUNIT_TEST(testSimplePolygons);
> + CPPUNIT_TEST(testComplexPolyhedrons);
> + CPPUNIT_TEST(testComplexPolygons);
> + CPPUNIT_TEST_SUITE_END();
> +
> +
> + void almost_equal_points(const Point & p1, const Point & p2, double delta)
> + {
> + CPPUNIT_ASSERT_DOUBLES_EQUAL(p1.x(),p2.x(),delta);
> + CPPUNIT_ASSERT_DOUBLES_EQUAL(p1.y(),p2.y(),delta);
> + CPPUNIT_ASSERT_DOUBLES_EQUAL(p1.z(),p2.z(),delta);
> + }
> +
> + //Helper function to create reference polyhedrons.
> + void add_test_polyhedron(const Point_3 & p1,
> + const Point_3 & p2,
> + const Point_3 & p3,
> + const Point_3 & p4,
> + PolyhedronList & polyhedrons
> + )
> + {
> + Polyhedron_3 P;
> + P.make_tetrahedron(p1,p2,p3,p4);
> + Nef_polyhedron_3 N(P);
> + polyhedrons.push_back(N);
> + }
> +
> + void add_test_polyhedron(const Point_3 & p1,
> + const Point_3 & p2,
> + const Point_3 & p3,
> + PolyhedronList & polyhedrons
> + )
> + {
> + Polyhedron_3 P;
> + P.make_triangle(p1,p2,p3);
> + Nef_polyhedron_3 N(P);
> + polyhedrons.push_back(N);
> + }
> +
> + //Helper function to union disjoint polyhedrons and to compute the volume and barycenter.
> + //Indices indicate which polyhedrons should be unioned. No checks at all (index, disjointness etc)
> + //Computed polyhedrons, volumes and barycenters will be append to the given list.
> + void add_disjoint_polyhedrons(const IntList indices,
> + PolyhedronList & polyhedrons,
> + DoubleList & volumes,
> + PointList & points)
> + {
> + double volume = 0;
> + Point point(0,0,0);
> + Nef_polyhedron_3 polyhedron;
> +
> + for (IntListIterator i = indices.begin(); i != indices.end(); ++i)
> + {
> + polyhedron += polyhedrons[*i];
> + point += volumes[*i] * points[*i] ;
> + volume += volumes[*i];
> + }
> + point /= volume;
> +
> + polyhedrons.push_back(polyhedron);
> + volumes.push_back(volume);
> + points.push_back(point);
> + }
> +
> +
> + void testSimplePolyhedrons()
> + {
> + //Create origin and unit vectors.
> + Point_3 e0(0.0, 0.0, 0.0);
> +
> + Point_3 e1(1.0, 0.0, 0.0);
> + Point_3 e2(0.0, 1.0, 0.0);
> + Point_3 e3(0.0, 0.0, 1.0);
> +
> + Point_3 _e1(-1.0, 0.0, 0.0);
> + Point_3 _e2(0.0, -1.0, 0.0);
> + Point_3 _e3(0.0, 0.0, -1.0);
> +
> + //Create tetrahedrons with unit vectors and reference results.
> + DoubleList reference_volumes;
> + PointList reference_bary_centers;
> + PolyhedronList reference_polyhedrons;
> +
> + add_test_polyhedron(e0,e1,e2,e3,reference_polyhedrons);
> + reference_volumes.push_back(1.0/6.0);
> + reference_bary_centers.push_back(Point(0.25,0.25,0.25));
> +
> + add_test_polyhedron(e0,_e1,e2,e3,reference_polyhedrons);
> + reference_volumes.push_back(1.0/6.0);
> + reference_bary_centers.push_back(Point(-0.25,0.25,0.25));
> +
> + add_test_polyhedron(e0,e1,_e2,e3,reference_polyhedrons);
> + reference_volumes.push_back(1.0/6.0);
> + reference_bary_centers.push_back(Point(0.25,-0.25,0.25));
> +
> + add_test_polyhedron(e0,_e1,_e2,e3,reference_polyhedrons);
> + reference_volumes.push_back(1.0/6.0);
> + reference_bary_centers.push_back(Point(-0.25,-0.25,0.25));
> +
> + add_test_polyhedron(e0,e1,e2,_e3,reference_polyhedrons);
> + reference_volumes.push_back(1.0/6.0);
> + reference_bary_centers.push_back(Point(0.25,0.25,-0.25));
> +
> + add_test_polyhedron(e0,_e1,e2,_e3,reference_polyhedrons);
> + reference_volumes.push_back(1.0/6.0);
> + reference_bary_centers.push_back(Point(-0.25,0.25,-0.25));
> +
> + add_test_polyhedron(e0,e1,_e2,_e3,reference_polyhedrons);
> + reference_volumes.push_back(1.0/6.0);
> + reference_bary_centers.push_back(Point(0.25,-0.25,-0.25));
> +
> + add_test_polyhedron(e0,_e1,_e2,_e3,reference_polyhedrons);
> + reference_volumes.push_back(1.0/6.0);
> + reference_bary_centers.push_back(Point(-0.25,-0.25,-0.25));
> +
> + //Add sum of polyhedrons
> + IntList add_indices;
> + add_indices.push_back(0);
> +
> + add_indices.push_back(1);
> + add_disjoint_polyhedrons(add_indices,
> + reference_polyhedrons,
> + reference_volumes,
> + reference_bary_centers);
> +
> + add_indices.push_back(2);
> + add_disjoint_polyhedrons(add_indices,
> + reference_polyhedrons,
> + reference_volumes,
> + reference_bary_centers);
> +
> + add_indices.push_back(3);
> + add_disjoint_polyhedrons(add_indices,
> + reference_polyhedrons,
> + reference_volumes,
> + reference_bary_centers);
> +
> + add_indices.push_back(4);
> + add_disjoint_polyhedrons(add_indices,
> + reference_polyhedrons,
> + reference_volumes,
> + reference_bary_centers);
> +
> + add_indices.push_back(5);
> + add_disjoint_polyhedrons(add_indices,
> + reference_polyhedrons,
> + reference_volumes,
> + reference_bary_centers);
> +
> + add_indices.push_back(6);
> + add_disjoint_polyhedrons(add_indices,
> + reference_polyhedrons,
> + reference_volumes,
> + reference_bary_centers);
> +
> + add_indices.push_back(7);
> + add_disjoint_polyhedrons(add_indices,
> + reference_polyhedrons,
> + reference_volumes,
> + reference_bary_centers);
> +
> + //Add translated version
> + //Upper halfspace
> + Nef_polyhedron_3 polyhedron = reference_polyhedrons[0];
> + polyhedron.transform(Aff_transformation_3(CGAL::TRANSLATION, Vector_3(1, 1, 1)));
> + reference_polyhedrons.push_back(polyhedron);
> + reference_volumes.push_back(1.0/6.0);
> + reference_bary_centers.push_back(Point(1.25,1.25,1.25));
> +
> + polyhedron = reference_polyhedrons[1];
> + polyhedron.transform(Aff_transformation_3(CGAL::TRANSLATION, Vector_3(-1, 1, 1)));
> + reference_polyhedrons.push_back(polyhedron);
> + reference_volumes.push_back(1.0/6.0);
> + reference_bary_centers.push_back(Point(-1.25,1.25,1.25));
> +
> + polyhedron = reference_polyhedrons[2];
> + polyhedron.transform(Aff_transformation_3(CGAL::TRANSLATION, Vector_3(1, -1, 1)));
> + reference_polyhedrons.push_back(polyhedron);
> + reference_volumes.push_back(1.0/6.0);
> + reference_bary_centers.push_back(Point(1.25,-1.25,1.25));
> +
> + polyhedron = reference_polyhedrons[3];
> + polyhedron.transform(Aff_transformation_3(CGAL::TRANSLATION, Vector_3(-1, -1, 1)));
> + reference_polyhedrons.push_back(polyhedron);
> + reference_volumes.push_back(1.0/6.0);
> + reference_bary_centers.push_back(Point(-1.25,-1.25,1.25));
> +
> + polyhedron = reference_polyhedrons[4];
> + polyhedron.transform(Aff_transformation_3(CGAL::TRANSLATION, Vector_3(1, 1, -1)));
> + reference_polyhedrons.push_back(polyhedron);
> + reference_volumes.push_back(1.0/6.0);
> + reference_bary_centers.push_back(Point(1.25,1.25,-1.25));
> +
> + polyhedron = reference_polyhedrons[5];
> + polyhedron.transform(Aff_transformation_3(CGAL::TRANSLATION, Vector_3(-1, 1, -1)));
> + reference_polyhedrons.push_back(polyhedron);
> + reference_volumes.push_back(1.0/6.0);
> + reference_bary_centers.push_back(Point(-1.25,1.25,-1.25));
> +
> + polyhedron = reference_polyhedrons[6];
> + polyhedron.transform(Aff_transformation_3(CGAL::TRANSLATION, Vector_3(1, -1, -1)));
> + reference_polyhedrons.push_back(polyhedron);
> + reference_volumes.push_back(1.0/6.0);
> + reference_bary_centers.push_back(Point(1.25,-1.25,-1.25));
> +
> + polyhedron = reference_polyhedrons[7];
> + polyhedron.transform(Aff_transformation_3(CGAL::TRANSLATION, Vector_3(-1, -1, -1)));
> + reference_polyhedrons.push_back(polyhedron);
> + reference_volumes.push_back(1.0/6.0);
> + reference_bary_centers.push_back(Point(-1.25,-1.25,-1.25));
> +
> + //Add the disjoint union of the translated polyhedrons.
> + add_indices.clear();
> + add_indices.push_back(15);
> +
> + add_indices.push_back(16);
> + add_disjoint_polyhedrons(add_indices,
> + reference_polyhedrons,
> + reference_volumes,
> + reference_bary_centers);
> +
> + add_indices.push_back(17);
> + add_disjoint_polyhedrons(add_indices,
> + reference_polyhedrons,
> + reference_volumes,
> + reference_bary_centers);
> +
> + add_indices.push_back(18);
> + add_disjoint_polyhedrons(add_indices,
> + reference_polyhedrons,
> + reference_volumes,
> + reference_bary_centers);
> +
> + add_indices.push_back(19);
> + add_disjoint_polyhedrons(add_indices,
> + reference_polyhedrons,
> + reference_volumes,
> + reference_bary_centers);
> +
> + add_indices.push_back(20);
> + add_disjoint_polyhedrons(add_indices,
> + reference_polyhedrons,
> + reference_volumes,
> + reference_bary_centers);
> +
> + add_indices.push_back(21);
> + add_disjoint_polyhedrons(add_indices,
> + reference_polyhedrons,
> + reference_volumes,
> + reference_bary_centers);
> +
> + add_indices.push_back(22);
> + add_disjoint_polyhedrons(add_indices,
> + reference_polyhedrons,
> + reference_volumes,
> + reference_bary_centers);
> +
> + //Instantiate quadrature rule
> + BaryCenterQuadrature quadrature_rule;
> +
> + //Check volume and barycenter for polyhedrons
> + for (dolfin::uint i = 0; i < reference_polyhedrons.size(); ++i)
> + {
> + quadrature_rule.compute_quadrature_rule(reference_polyhedrons[i]);
> + CPPUNIT_ASSERT_DOUBLES_EQUAL(reference_volumes[i],
> + quadrature_rule.weights()[0], 1.0e-12);
> + almost_equal_points(reference_bary_centers[i],
> + quadrature_rule.points()[0], 1.0e-12);
> + }
> + }
> +
> + void testSimplePolygons()
> + {
> + //Create origin and unit vectors.
> + Point_3 e0(0.0, 0.0, 0.0);
> +
> + Point_3 e1(1.0, 0.0, 0.0);
> + Point_3 e2(0.0, 1.0, 0.0);
> + Point_3 e3(0.0, 0.0, 1.0);
> +
> + Point_3 _e1(-1.0, 0.0, 0.0);
> + Point_3 _e2(0.0, -1.0, 0.0);
> + Point_3 _e3(0.0, 0.0, -1.0);
> +
> + //Create tetrahedrons with unit vectors and reference results.
> + DoubleList reference_volumes;
> + PointList reference_bary_centers;
> + PolyhedronList reference_polyhedrons;
> +
> + //skew plane upper e3 plane
> + add_test_polyhedron(e1,e2,e3,reference_polyhedrons);
> + //todo find exact values
> + reference_volumes.push_back(8.660254e-01);
> + reference_bary_centers.push_back(Point(1.0/3.0,1.0/3.0,1.0/3.0));
> +
> + add_test_polyhedron(_e1,e2,e3,reference_polyhedrons);
> + //todo find exact values
> + reference_volumes.push_back(8.660254e-01);
> + reference_bary_centers.push_back(Point(-1.0/3.0,1.0/3.0,1.0/3.0));
> +
> + add_test_polyhedron(_e1,_e2,e3,reference_polyhedrons);
> + //todo find exact values
> + reference_volumes.push_back(8.660254e-01);
> + reference_bary_centers.push_back(Point(-1.0/3.0,-1.0/3.0,1.0/3.0));
> +
> + add_test_polyhedron(e1,_e2,e3,reference_polyhedrons);
> + //todo find exact values
> + reference_volumes.push_back(8.660254e-01);
> + reference_bary_centers.push_back(Point(1.0/3.0,-1.0/3.0,1.0/3.0));
> +
> + //skew plane lower -e3 plane
> + add_test_polyhedron(e1,e2,_e3,reference_polyhedrons);
> + //todo find exact values
> + reference_volumes.push_back(8.660254e-01);
> + reference_bary_centers.push_back(Point(1.0/3.0,1.0/3.0,-1.0/3.0));
> +
> + add_test_polyhedron(_e1,e2,_e3,reference_polyhedrons);
> + //todo find exact values
> + reference_volumes.push_back(8.660254e-01);
> + reference_bary_centers.push_back(Point(-1.0/3.0,1.0/3.0,-1.0/3.0));
> +
> + add_test_polyhedron(_e1,_e2,_e3,reference_polyhedrons);
> + //todo find exact values
> + reference_volumes.push_back(8.660254e-01);
> + reference_bary_centers.push_back(Point(-1.0/3.0,-1.0/3.0,-1.0/3.0));
> +
> + add_test_polyhedron(e1,_e2,_e3,reference_polyhedrons);
> + //todo find exact values
> + reference_volumes.push_back(8.660254e-01);
> + reference_bary_centers.push_back(Point(1.0/3.0,-1.0/3.0,-1.0/3.0));
> +
> + //e1-e2 plane
> + add_test_polyhedron(e0,e1,e2,reference_polyhedrons);
> + reference_volumes.push_back(1.0/2.0);
> + reference_bary_centers.push_back(Point(1.0/3.0,1.0/3.0,0.0));
> +
> + add_test_polyhedron(e0,_e1,e2,reference_polyhedrons);
> + reference_volumes.push_back(1.0/2.0);
> + reference_bary_centers.push_back(Point(-1.0/3.0,1.0/3.0,0.0));
> +
> + add_test_polyhedron(e0,_e1,_e2,reference_polyhedrons);
> + reference_volumes.push_back(1.0/2.0);
> + reference_bary_centers.push_back(Point(-1.0/3.0,-1.0/3.0,0.0));
> +
> + add_test_polyhedron(e0,e1,_e2,reference_polyhedrons);
> + reference_volumes.push_back(1.0/2.0);
> + reference_bary_centers.push_back(Point(1.0/3.0,-1.0/3.0,0.0));
> +
> + //e1-e3 plane
> + add_test_polyhedron(e0,e1,e3,reference_polyhedrons);
> + reference_volumes.push_back(1.0/2.0);
> + reference_bary_centers.push_back(Point(1.0/3.0,0.0,1.0/3.0));
> +
> + add_test_polyhedron(e0,_e1,e3,reference_polyhedrons);
> + reference_volumes.push_back(1.0/2.0);
> + reference_bary_centers.push_back(Point(-1.0/3.0,0.0,1.0/3.0));
> +
> + add_test_polyhedron(e0,_e1,_e3,reference_polyhedrons);
> + reference_volumes.push_back(1.0/2.0);
> + reference_bary_centers.push_back(Point(-1.0/3.0,0.0,-1.0/3.0));
> +
> + add_test_polyhedron(e0,e1,_e3,reference_polyhedrons);
> + reference_volumes.push_back(1.0/2.0);
> + reference_bary_centers.push_back(Point(1.0/3.0,0.0,-1.0/3.0));
> +
> + //e2-e3 plane
> + add_test_polyhedron(e0,e2,e3,reference_polyhedrons);
> + reference_volumes.push_back(1.0/2.0);
> + reference_bary_centers.push_back(Point(0.0,1.0/3.0,1.0/3.0));
> +
> + add_test_polyhedron(e0,_e2,e3,reference_polyhedrons);
> + reference_volumes.push_back(1.0/2.0);
> + reference_bary_centers.push_back(Point(0.0,-1.0/3.0,1.0/3.0));
> +
> + add_test_polyhedron(e0,_e2,_e3,reference_polyhedrons);
> + reference_volumes.push_back(1.0/2.0);
> + reference_bary_centers.push_back(Point(0.0,-1.0/3.0,-1.0/3.0));
> +
> + add_test_polyhedron(e0,e2,_e3,reference_polyhedrons);
> + reference_volumes.push_back(1.0/2.0);
> + reference_bary_centers.push_back(Point(0.0,1.0/3.0,-1.0/3.0));
> +
> + //Test sum of polyhedrons
> + IntList add_indices;
> + add_indices.push_back(0);
> +
> +
> + add_indices.push_back(1);
> + add_disjoint_polyhedrons(add_indices,
> + reference_polyhedrons,
> + reference_volumes,
> + reference_bary_centers);
> +
> + add_indices.push_back(2);
> + add_disjoint_polyhedrons(add_indices,
> + reference_polyhedrons,
> + reference_volumes,
> + reference_bary_centers);
> +
> + add_indices.push_back(3);
> + add_disjoint_polyhedrons(add_indices,
> + reference_polyhedrons,
> + reference_volumes,
> + reference_bary_centers);
> +
> + add_indices.push_back(4);
> + add_disjoint_polyhedrons(add_indices,
> + reference_polyhedrons,
> + reference_volumes,
> + reference_bary_centers);
> +
> + add_indices.push_back(5);
> + add_disjoint_polyhedrons(add_indices,
> + reference_polyhedrons,
> + reference_volumes,
> + reference_bary_centers);
> +
> + add_indices.push_back(6);
> + add_disjoint_polyhedrons(add_indices,
> + reference_polyhedrons,
> + reference_volumes,
> + reference_bary_centers);
> +
> + add_indices.push_back(7);
> + add_disjoint_polyhedrons(add_indices,
> + reference_polyhedrons,
> + reference_volumes,
> + reference_bary_centers);
> +
> + add_indices.push_back(8);
> + add_disjoint_polyhedrons(add_indices,
> + reference_polyhedrons,
> + reference_volumes,
> + reference_bary_centers);
> +
> + add_indices.push_back(9);
> + add_disjoint_polyhedrons(add_indices,
> + reference_polyhedrons,
> + reference_volumes,
> + reference_bary_centers);
> +
> +
> + add_indices.push_back(10);
> + add_disjoint_polyhedrons(add_indices,
> + reference_polyhedrons,
> + reference_volumes,
> + reference_bary_centers);
> +
> + add_indices.push_back(11);
> + add_disjoint_polyhedrons(add_indices,
> + reference_polyhedrons,
> + reference_volumes,
> + reference_bary_centers);
> +
> + add_indices.push_back(12);
> + add_disjoint_polyhedrons(add_indices,
> + reference_polyhedrons,
> + reference_volumes,
> + reference_bary_centers);
> +
> + add_indices.push_back(13);
> + add_disjoint_polyhedrons(add_indices,
> + reference_polyhedrons,
> + reference_volumes,
> + reference_bary_centers);
> +
> + add_indices.push_back(14);
> + add_disjoint_polyhedrons(add_indices,
> + reference_polyhedrons,
> + reference_volumes,
> + reference_bary_centers);
> +
> + add_indices.push_back(15);
> + add_disjoint_polyhedrons(add_indices,
> + reference_polyhedrons,
> + reference_volumes,
> + reference_bary_centers);
> +
> + Nef_polyhedron_3 polyhedron = reference_polyhedrons[0];
> + polyhedron.transform(Aff_transformation_3(CGAL::TRANSLATION, Vector_3(1, 1, 1)));
> + reference_polyhedrons.push_back(polyhedron);
> + reference_volumes.push_back(reference_volumes[0]);
> + reference_bary_centers.push_back(reference_bary_centers[0] + Point(1,1,1));
> +
> + polyhedron = reference_polyhedrons[1];
> + polyhedron.transform(Aff_transformation_3(CGAL::TRANSLATION, Vector_3(-1, 1, 1)));
> + reference_polyhedrons.push_back(polyhedron);
> + reference_volumes.push_back(reference_volumes[1]);
> + reference_bary_centers.push_back(reference_bary_centers[1] + Point(-1,1,1));
> +
> + polyhedron = reference_polyhedrons[2];
> + polyhedron.transform(Aff_transformation_3(CGAL::TRANSLATION, Vector_3(1, -1, 1)));
> + reference_polyhedrons.push_back(polyhedron);
> + reference_volumes.push_back(reference_volumes[2]);
> + reference_bary_centers.push_back(reference_bary_centers[2] + Point(1,-1,1));
> +
> + polyhedron = reference_polyhedrons[3];
> + polyhedron.transform(Aff_transformation_3(CGAL::TRANSLATION, Vector_3(-1, -1, 1)));
> + reference_polyhedrons.push_back(polyhedron);
> + reference_volumes.push_back(reference_volumes[3]);
> + reference_bary_centers.push_back(reference_bary_centers[3] + Point(-1,-1,1));
> +
> + polyhedron = reference_polyhedrons[4];
> + polyhedron.transform(Aff_transformation_3(CGAL::TRANSLATION, Vector_3(1, 1, -1)));
> + reference_polyhedrons.push_back(polyhedron);
> + reference_volumes.push_back(reference_volumes[4]);
> + reference_bary_centers.push_back(reference_bary_centers[4] + Point(1,1,-1));
> +
> + polyhedron = reference_polyhedrons[5];
> + polyhedron.transform(Aff_transformation_3(CGAL::TRANSLATION, Vector_3(-1, 1, -1)));
> + reference_polyhedrons.push_back(polyhedron);
> + reference_volumes.push_back(reference_volumes[5]);
> + reference_bary_centers.push_back(reference_bary_centers[5] + Point(-1,1,-1));
> +
> + polyhedron = reference_polyhedrons[6];
> + polyhedron.transform(Aff_transformation_3(CGAL::TRANSLATION, Vector_3(1, -1, -1)));
> + reference_polyhedrons.push_back(polyhedron);
> + reference_volumes.push_back(reference_volumes[6]);
> + reference_bary_centers.push_back(reference_bary_centers[6] + Point(1,-1,-1));
> +
> + polyhedron = reference_polyhedrons[7];
> + polyhedron.transform(Aff_transformation_3(CGAL::TRANSLATION, Vector_3(-1, -1, -1)));
> + reference_polyhedrons.push_back(polyhedron);
> + reference_volumes.push_back(reference_volumes[7]);
> + reference_bary_centers.push_back(reference_bary_centers[7] + Point(-1,-1,-1));
> +
> + //Instantiate quadrature rule
> + BaryCenterQuadrature quadrature_rule;
> +
> + //Check volume and barycenter for polyhedrons
> + for (dolfin::uint i = 0; i < reference_polyhedrons.size(); ++i)
> + {
> + quadrature_rule.compute_quadrature_rule(reference_polyhedrons[i]);
> + CPPUNIT_ASSERT_DOUBLES_EQUAL(reference_volumes[i],
> + quadrature_rule.weights()[0], 1.0e-5);
> + almost_equal_points(reference_bary_centers[i],
> + quadrature_rule.points()[0], 1.0e-5);
> + }
> + }
> +
> + void testComplexPolyhedrons()
> + {
> + }
> +
> + void testComplexPolygons()
> + {
> + }
> +
> +};
> +
> +CPPUNIT_TEST_SUITE_REGISTRATION(BaryCenter);
> +
> +int main()
> +{
> + DOLFIN_TEST;
> +}
>
> === added file 'nitsche_method/integration_on_complex_domains/gauss_projection_green/tetra'
> --- nitsche_method/integration_on_complex_domains/gauss_projection_green/tetra 1970-01-01 00:00:00 +0000
> +++ nitsche_method/integration_on_complex_domains/gauss_projection_green/tetra 2010-03-18 18:15:26 +0000
> @@ -0,0 +1,16 @@
> +4
> +
> +0 0 0
> +5 0 0
> +0 4 0
> +0 0 3
> +
> +4
> +
> +3 0 3 2
> +3 3 0 1
> +3 2 1 0
> +3 1 2 3
> +
> +
> +
>
> === added file 'nitsche_method/integration_on_complex_domains/gauss_projection_green/volInt.c'
> --- nitsche_method/integration_on_complex_domains/gauss_projection_green/volInt.c 1970-01-01 00:00:00 +0000
> +++ nitsche_method/integration_on_complex_domains/gauss_projection_green/volInt.c 2010-03-26 13:05:51 +0000
> @@ -0,0 +1,383 @@
> +
> +
> + /*******************************************************
> + * *
> + * volInt.c *
> + * *
> + * This code computes volume integrals needed for *
> + * determining mass properties of polyhedral bodies. *
> + * *
> + * For more information, see the accompanying README *
> + * file, and the paper *
> + * *
> + * Brian Mirtich, "Fast and Accurate Computation of *
> + * Polyhedral Mass Properties," journal of graphics *
> + * tools, volume 1, number 1, 1996. *
> + * *
> + * This source code is public domain, and may be used *
> + * in any way, shape or form, free of charge. *
> + * *
> + * Copyright 1995 by Brian Mirtich *
> + * *
> + * mirtich@xxxxxxxxxxxxxxx *
> + * http://www.cs.berkeley.edu/~mirtich *
> + * *
> + *******************************************************/
> +
> +/*
> + Revision history
> +
> + 26 Jan 1996 Program creation.
> +
> + 3 Aug 1996 Corrected bug arising when polyhedron density
> + is not 1.0. Changes confined to function main().
> + Thanks to Zoran Popovic for catching this one.
> +
> + 27 May 1997 Corrected sign error in translation of inertia
> + product terms to center of mass frame. Changes
> + confined to function main(). Thanks to
> + Chris Hecker.
> +*/
> +
> +
> +
> +#include <stdio.h>
> +#include <string.h>
> +#include <math.h>
> +
> +using namespace std;
> +
> +/*
> + ============================================================================
> + constants
> + ============================================================================
> +*/
> +
> +#define MAX_VERTS 100 /* maximum number of polyhedral vertices */
> +#define MAX_FACES 100 /* maximum number of polyhedral faces */
> +#define MAX_POLYGON_SZ 10 /* maximum number of verts per polygonal face */
> +
> +#define X 0
> +#define Y 1
> +#define Z 2
> +
> +/*
> + ============================================================================
> + macros
> + ============================================================================
> +*/
> +
> +#define SQR(x) ((x)*(x))
> +#define CUBE(x) ((x)*(x)*(x))
> +
> +/*
> + ============================================================================
> + data structures
> + ============================================================================
> +*/
> +
> +typedef struct {
> + int numVerts;
> + double norm[3];
> + double w;
> + int verts[MAX_POLYGON_SZ];
> + struct polyhedron *poly;
> +} FACE;
> +
> +typedef struct polyhedron {
> + int numVerts, numFaces;
> + double verts[MAX_VERTS][3];
> + FACE faces[MAX_FACES];
> +} POLYHEDRON;
> +
> +
> +/*
> + ============================================================================
> + globals
> + ============================================================================
> +*/
> +
> +static int A; /* alpha */
> +static int B; /* beta */
> +static int C; /* gamma */
> +
> +/* projection integrals */
> +static double P1, Pa, Pb, Paa, Pab, Pbb, Paaa, Paab, Pabb, Pbbb;
> +
> +/* face integrals */
> +static double Fa, Fb, Fc, Faa, Fbb, Fcc, Faaa, Fbbb, Fccc, Faab, Fbbc, Fcca;
> +
> +/* volume integrals */
> +static double T0, T1[3], T2[3], TP[3];
> +
> +
> +/*
> + ============================================================================
> + read in a polyhedron
> + ============================================================================
> +*/
> +
> +void readPolyhedron(char *name, POLYHEDRON *p)
> +{
> + FILE *fp;
> + char line[200], *c;
> + int i, j, n;
> + double dx1, dy1, dz1, dx2, dy2, dz2, nx, ny, nz, len;
> + FACE *f;
> +
> +
> +/* if (!(fp = fopen(name, "r"))) {*/
> +/* printf("i/o error\n");*/
> +/* exit(1);*/
> + }
> +
> + fscanf(fp, "%d", &p->numVerts);
> + printf("Reading in %d vertices\n", p->numVerts);
> + for (i = 0; i < p->numVerts; i++)
> + fscanf(fp, "%lf %lf %lf",
> + &p->verts[i][X], &p->verts[i][Y], &p->verts[i][Z]);
> +
> + fscanf(fp, "%d", &p->numFaces);
> + printf("Reading in %d faces\n", p->numFaces);
> + for (i = 0; i < p->numFaces; i++) {
> + f = &p->faces[i];
> + f->poly = p;
> + fscanf(fp, "%d", &f->numVerts);
> + for (j = 0; j < f->numVerts; j++) fscanf(fp, "%d", &f->verts[j]);
> +
> + /* compute face normal and offset w from first 3 vertices */
> + dx1 = p->verts[f->verts[1]][X] - p->verts[f->verts[0]][X];
> + dy1 = p->verts[f->verts[1]][Y] - p->verts[f->verts[0]][Y];
> + dz1 = p->verts[f->verts[1]][Z] - p->verts[f->verts[0]][Z];
> + dx2 = p->verts[f->verts[2]][X] - p->verts[f->verts[1]][X];
> + dy2 = p->verts[f->verts[2]][Y] - p->verts[f->verts[1]][Y];
> + dz2 = p->verts[f->verts[2]][Z] - p->verts[f->verts[1]][Z];
> + nx = dy1 * dz2 - dy2 * dz1;
> + ny = dz1 * dx2 - dz2 * dx1;
> + nz = dx1 * dy2 - dx2 * dy1;
> + len = sqrt(nx * nx + ny * ny + nz * nz);
> + f->norm[X] = nx / len;
> + f->norm[Y] = ny / len;
> + f->norm[Z] = nz / len;
> + f->w = - f->norm[X] * p->verts[f->verts[0]][X]
> + - f->norm[Y] * p->verts[f->verts[0]][Y]
> + - f->norm[Z] * p->verts[f->verts[0]][Z];
> +
> + }
> +
> + fclose(fp);
> +
> +}
> +
> +/*
> + ============================================================================
> + compute mass properties
> + ============================================================================
> +*/
> +
> +
> +/* compute various integrations over projection of face */
> +void compProjectionIntegrals(FACE *f)
> +{
> + double a0, a1, da;
> + double b0, b1, db;
> + double a0_2, a0_3, a0_4, b0_2, b0_3, b0_4;
> + double a1_2, a1_3, b1_2, b1_3;
> + double C1, Ca, Caa, Caaa, Cb, Cbb, Cbbb;
> + double Cab, Kab, Caab, Kaab, Cabb, Kabb;
> + int i;
> +
> + P1 = Pa = Pb = Paa = Pab = Pbb = Paaa = Paab = Pabb = Pbbb = 0.0;
> +
> + for (i = 0; i < f->numVerts; i++) {
> + a0 = f->poly->verts[f->verts[i]][A];
> + b0 = f->poly->verts[f->verts[i]][B];
> + a1 = f->poly->verts[f->verts[(i+1) % f->numVerts]][A];
> + b1 = f->poly->verts[f->verts[(i+1) % f->numVerts]][B];
> + da = a1 - a0;
> + db = b1 - b0;
> + a0_2 = a0 * a0; a0_3 = a0_2 * a0; a0_4 = a0_3 * a0;
> + b0_2 = b0 * b0; b0_3 = b0_2 * b0; b0_4 = b0_3 * b0;
> + a1_2 = a1 * a1; a1_3 = a1_2 * a1;
> + b1_2 = b1 * b1; b1_3 = b1_2 * b1;
> +
> + C1 = a1 + a0;
> + Ca = a1*C1 + a0_2; Caa = a1*Ca + a0_3; Caaa = a1*Caa + a0_4;
> + Cb = b1*(b1 + b0) + b0_2; Cbb = b1*Cb + b0_3; Cbbb = b1*Cbb + b0_4;
> + Cab = 3*a1_2 + 2*a1*a0 + a0_2; Kab = a1_2 + 2*a1*a0 + 3*a0_2;
> + Caab = a0*Cab + 4*a1_3; Kaab = a1*Kab + 4*a0_3;
> + Cabb = 4*b1_3 + 3*b1_2*b0 + 2*b1*b0_2 + b0_3;
> + Kabb = b1_3 + 2*b1_2*b0 + 3*b1*b0_2 + 4*b0_3;
> +
> + P1 += db*C1;
> + Pa += db*Ca;
> + Paa += db*Caa;
> + Paaa += db*Caaa;
> + Pb += da*Cb;
> + Pbb += da*Cbb;
> + Pbbb += da*Cbbb;
> + Pab += db*(b1*Cab + b0*Kab);
> + Paab += db*(b1*Caab + b0*Kaab);
> + Pabb += da*(a1*Cabb + a0*Kabb);
> + }
> +
> + P1 /= 2.0;
> + Pa /= 6.0;
> + Paa /= 12.0;
> + Paaa /= 20.0;
> + Pb /= -6.0;
> + Pbb /= -12.0;
> + Pbbb /= -20.0;
> + Pab /= 24.0;
> + Paab /= 60.0;
> + Pabb /= -60.0;
> +}
> +
> +compFaceIntegrals(FACE *f)
> +{
> + double *n, w;
> + double k1, k2, k3, k4;
> +
> + compProjectionIntegrals(f);
> +
> + w = f->w;
> + n = f->norm;
> + k1 = 1 / n[C]; k2 = k1 * k1; k3 = k2 * k1; k4 = k3 * k1;
> +
> + Fa = k1 * Pa;
> + Fb = k1 * Pb;
> + Fc = -k2 * (n[A]*Pa + n[B]*Pb + w*P1);
> +
> + Faa = k1 * Paa;
> + Fbb = k1 * Pbb;
> + Fcc = k3 * (SQR(n[A])*Paa + 2*n[A]*n[B]*Pab + SQR(n[B])*Pbb
> + + w*(2*(n[A]*Pa + n[B]*Pb) + w*P1));
> +
> + Faaa = k1 * Paaa;
> + Fbbb = k1 * Pbbb;
> + Fccc = -k4 * (CUBE(n[A])*Paaa + 3*SQR(n[A])*n[B]*Paab
> + + 3*n[A]*SQR(n[B])*Pabb + CUBE(n[B])*Pbbb
> + + 3*w*(SQR(n[A])*Paa + 2*n[A]*n[B]*Pab + SQR(n[B])*Pbb)
> + + w*w*(3*(n[A]*Pa + n[B]*Pb) + w*P1));
> +
> + Faab = k1 * Paab;
> + Fbbc = -k2 * (n[A]*Pabb + n[B]*Pbbb + w*Pbb);
> + Fcca = k3 * (SQR(n[A])*Paaa + 2*n[A]*n[B]*Paab + SQR(n[B])*Pabb
> + + w*(2*(n[A]*Paa + n[B]*Pab) + w*Pa));
> +}
> +
> +void compVolumeIntegrals(POLYHEDRON *p)
> +{
> + FACE *f;
> + double nx, ny, nz;
> + int i;
> +
> + T0 = T1[X] = T1[Y] = T1[Z]
> + = T2[X] = T2[Y] = T2[Z]
> + = TP[X] = TP[Y] = TP[Z] = 0;
> +
> + for (i = 0; i < p->numFaces; i++) {
> +
> + f = &p->faces[i];
> +
> + nx = fabs(f->norm[X]);
> + ny = fabs(f->norm[Y]);
> + nz = fabs(f->norm[Z]);
> + if (nx > ny && nx > nz) C = X;
> + else C = (ny > nz) ? Y : Z;
> + A = (C + 1) % 3;
> + B = (A + 1) % 3;
> +
> + compFaceIntegrals(f);
> +
> + T0 += f->norm[X] * ((A == X) ? Fa : ((B == X) ? Fb : Fc));
> +
> + T1[A] += f->norm[A] * Faa;
> + T1[B] += f->norm[B] * Fbb;
> + T1[C] += f->norm[C] * Fcc;
> + T2[A] += f->norm[A] * Faaa;
> + T2[B] += f->norm[B] * Fbbb;
> + T2[C] += f->norm[C] * Fccc;
> + TP[A] += f->norm[A] * Faab;
> + TP[B] += f->norm[B] * Fbbc;
> + TP[C] += f->norm[C] * Fcca;
> + }
> +
> + T1[X] /= 2; T1[Y] /= 2; T1[Z] /= 2;
> + T2[X] /= 3; T2[Y] /= 3; T2[Z] /= 3;
> + TP[X] /= 2; TP[Y] /= 2; TP[Z] /= 2;
> +}
> +
> +
> +/*
> + ============================================================================
> + main
> + ============================================================================
> +*/
> +
> +
> +int main(int argc, char *argv[])
> +{
> + POLYHEDRON p;
> + double density, mass;
> + double r[3]; /* center of mass */
> + double J[3][3]; /* inertia tensor */
> +
> + if (argc != 2) {
> + printf("usage: %s <polyhedron geometry filename>\n", argv[0]);
> +/* exit(0);*/
> + return 0;
> + }
> +
> + readPolyhedron(argv[1], &p);
> +
> + compVolumeIntegrals(&p);
> +
> +
> + printf("\nT1 = %+20.6f\n\n", T0);
> +
> + printf("Tx = %+20.6f\n", T1[X]);
> + printf("Ty = %+20.6f\n", T1[Y]);
> + printf("Tz = %+20.6f\n\n", T1[Z]);
> +
> + printf("Txx = %+20.6f\n", T2[X]);
> + printf("Tyy = %+20.6f\n", T2[Y]);
> + printf("Tzz = %+20.6f\n\n", T2[Z]);
> +
> + printf("Txy = %+20.6f\n", TP[X]);
> + printf("Tyz = %+20.6f\n", TP[Y]);
> + printf("Tzx = %+20.6f\n\n", TP[Z]);
> +
> + density = 1.0; /* assume unit density */
> +
> + mass = density * T0;
> +
> + /* compute center of mass */
> + r[X] = T1[X] / T0;
> + r[Y] = T1[Y] / T0;
> + r[Z] = T1[Z] / T0;
> +
> + /* compute inertia tensor */
> + J[X][X] = density * (T2[Y] + T2[Z]);
> + J[Y][Y] = density * (T2[Z] + T2[X]);
> + J[Z][Z] = density * (T2[X] + T2[Y]);
> + J[X][Y] = J[Y][X] = - density * TP[X];
> + J[Y][Z] = J[Z][Y] = - density * TP[Y];
> + J[Z][X] = J[X][Z] = - density * TP[Z];
> +
> + /* translate inertia tensor to center of mass */
> + J[X][X] -= mass * (r[Y]*r[Y] + r[Z]*r[Z]);
> + J[Y][Y] -= mass * (r[Z]*r[Z] + r[X]*r[X]);
> + J[Z][Z] -= mass * (r[X]*r[X] + r[Y]*r[Y]);
> + J[X][Y] = J[Y][X] += mass * r[X] * r[Y];
> + J[Y][Z] = J[Z][Y] += mass * r[Y] * r[Z];
> + J[Z][X] = J[X][Z] += mass * r[Z] * r[X];
> +
> + printf("center of mass: (%+12.6f,%+12.6f,%+12.6f)\n\n", r[X], r[Y], r[Z]);
> +
> + printf("inertia tensor with origin at c.o.m. :\n");
> + printf("%+15.6f %+15.6f %+15.6f\n", J[X][X], J[X][Y], J[X][Z]);
> + printf("%+15.6f %+15.6f %+15.6f\n", J[Y][X], J[Y][Y], J[Y][Z]);
> + printf("%+15.6f %+15.6f %+15.6f\n\n", J[Z][X], J[Z][Y], J[Z][Z]);
> +
> +}
>
> === added directory 'nitsche_method/integration_on_complex_domains/monte_carlo'
> === added directory 'nitsche_method/integration_on_complex_domains/triangulation'
> === added directory 'nitsche_method/overlapping_meshes'
> === added file 'nitsche_method/overlapping_meshes/MeshIntersector.cpp'
> --- nitsche_method/overlapping_meshes/MeshIntersector.cpp 1970-01-01 00:00:00 +0000
> +++ nitsche_method/overlapping_meshes/MeshIntersector.cpp 2010-03-16 14:39:09 +0000
> @@ -0,0 +1,175 @@
> +// =====================================================================================
> +//
> +// Copyright (C) 2010-01-16 André Massing
> +// Licensed under the GNU LGPL Version 2.1.
> +//
> +// Modified by André Massing, 2010
> +//
> +// First added: 2010-01-16
> +// back changed: 2010-01-27
> +//
> +//Author: André Massing (am), massing@xxxxxxxxx
> +//Company: Simula Research Laboratory, Fornebu, Norway
> +//
> +// =====================================================================================
> +
> +#include "OverlappingMeshes.h"
> +
> +#include <dolfin/mesh/BoundaryMesh.h>
> +#include <dolfin/mesh/Vertex.h>
> +#include <dolfin/mesh/Cell.h>
> +
> +#include <dolfin/log/dolfin_log.h>
> +#include <dolfin/common/NoDeleter.h>
> +
> +#include <iostream>
> +
> +using namespace dolfin;
> +using dolfin::uint;
> +
> +
> + MeshIntersection::OverlapData::OverlapData(boost::shared_ptr<const Mesh> mesh)
> +: mesh(mesh)
> +{
> + intersected_domain = MeshFunction<uint>(*mesh,mesh->topology().dim());
> + intersected_domain = 0;
> +}
> +
> +MeshIntersection::MeshIntersection(const Mesh & mesh_1, const Mesh & mesh_2)
> + :
> + _overlapped_domain(new OverlapData(reference_to_no_delete_pointer(mesh_1)))
> + ,_overlapping_domain(new OverlapData(reference_to_no_delete_pointer(mesh_2)))
> + ,_overlapped_boundary(new OverlapData(boost::shared_ptr<const Mesh>(new BoundaryMesh(mesh_1))))
> + ,_overlapping_boundary(new OverlapData(boost::shared_ptr<const Mesh>(new BoundaryMesh(mesh_2))))
> +{}
> +
> +//MeshIntersection::MeshIntersection(boost::shared_ptr<const Mesh> mesh_1, boost::shared_ptr<const Mesh> mesh_2)
> +// :
> +// _overlapped_domain(mesh_1)
> +// ,_overlapping_domain(mesh_2)
> +// _overlapped_boundary(boost::shared_ptr<const Mesh>(new BoundaryMesh(mesh_1)))
> +// ,_overlapping_boundary(boost::shared_ptr<const Mesh>(new BoundaryMesh(mesh_2)))
> +//{}
> +
> +
> +void MeshIntersection::compute_overlap_map()
> +{
> +
> + const Mesh & mesh_1 = *(_overlapped_domain->mesh);
> + const Mesh & mesh_2 = *(_overlapping_domain->mesh);
> +
> + const Mesh & boundary_1 = *(_overlapped_boundary->mesh);
> + const Mesh & boundary_2 = *(_overlapping_boundary->mesh);
> +
> + EntityEntitiesMap & cell_cell_map = _overlapped_domain->entity_entities_map;
> + EntityEntitiesMap & facet_1_cell_2_map = _overlapped_boundary->entity_entities_map;
> + EntityEntitiesMap & facet_2_cell_1_map = _overlapping_boundary->entity_entities_map;
> +
> + CellIterator cut_cell(mesh_1);
> + CellIterator cut_facet(boundary_2);
> +
> + //Step 1:
> + //Intersect boundary of mesh_2 with mesh_1.
> + //to get the *partially*
> + //intersected cells of mesh_1.
> + //This calculates:
> + // a) *partially* intersected cells of mesh1
> + // b) the exterior facets of mesh_2 which are (part) of the artificial interface.
> + // c) *partially* intersected exterior facets of mesh1 and mesh2.
> +
> + for (CellIterator cell(boundary_2); !cell.end(); ++cell)
> + {
> + mesh_1.all_intersected_entities(*cell, facet_2_cell_1_map[cell->index()]);
> +
> + if (facet_2_cell_1_map[cell->index()].empty())
> + facet_2_cell_1_map.erase(cell->index());
> + //If not empty add cell index and intersecting cell index to the map.
> + else
> + {
> + //Iterate of intersected cell of mesh1, find the overlapping cells of
> + //mesh 2 and mark cells in mesh1 as partially overlapped.
> + for (EntityListIter cell_iter = facet_2_cell_1_map[cell->index()].begin(); cell_iter != facet_2_cell_1_map[cell->index()].end(); ++cell_iter)
> + {
> + mesh_2.all_intersected_entities(cut_cell[*cell_iter], cell_cell_map[*cell_iter]);
> + _overlapped_domain->intersected_domain[*cell_iter] = 1;
> + }
> +
> + //Compute partially overlapped boundary cells of mesh1 and mesh2.
> + //
> + //@remark: Clarify whether it is faster to check first if any and then
> + //if compute indeces or just try to compute immediately and erase if the cell
> + //index container is empty. Rational: We want to avoid a copy operator
> + //(linear). A "any intersection" test should have the same complexity as
> + //a "compute all intersection" if no intersection occurs. If a
> + //intersection occurrs, we want to compute all intersection anyway.
> + //1. Version Compute right away and deleting (delete operation is amortized constant) if empty
> + //2. Version Check first and compute if intersected.
> + //3. Compute and assign if not empty
> + //@remark What is faster? Recompute intersecting cells from mesh2 for the
> + //exterior facets in mesh1 or map facet index to cell index, and assign
> + //their cell set (which might be bigger as the set, which really only
> + //intersects the facet).
> +
> + EntityList cut_faces;
> + boundary_1.all_intersected_entities(*cell,cut_faces);
> +
> + if (!cut_faces.empty())
> + {
> + _overlapping_boundary->intersected_domain[cell->index()] = 1;
> +
> + //Compute for each cut exterior facet in mesh1 the cutting cells in
> + //mesh2, mark facet as partially overlapped.
> + for (EntityListIter face_iter = cut_faces.begin(); face_iter != cut_faces.end(); ++face_iter)
> + {
> + mesh_2.all_intersected_entities(cut_facet[*face_iter],facet_1_cell_2_map[*face_iter]);
> + _overlapped_boundary->intersected_domain[*face_iter] = 1;
> + }
> + }
> + else
> + {
> + _overlapping_boundary->intersected_domain[cell->index()] = 2;
> + }
> + }
> + }
> +
> + //Step 2:
> + //Determine all cells of Mesh 1, which are fully overlapped. This is done by
> + //going through the cells, check if they are not partially overlapped and
> + //therefore must then be fully overlapped if any vertex is intersecting
> + //mesh2.
> + for (CellIterator cell(mesh_1); !cell.end(); ++cell)
> + {
> + if (_overlapped_domain->intersected_domain[cell->index()] != 1 && mesh_2.any_intersected_entity(VertexIterator(*cell)->point()) != -1)
> + _overlapped_domain->intersected_domain[cell->index()] = 2;
> + }
> +
> + //Step 3:
> + //Determine all cells of the boundary of mesh 1, which are fully overlapped.
> + //Same method as in Step 2.
> + for (CellIterator cell(boundary_1); !cell.end(); ++cell)
> + {
> + if (_overlapped_boundary->intersected_domain[cell->index()] != 1
> + && mesh_2.any_intersected_entity(VertexIterator(*cell)->point()) != -1)
> + _overlapped_boundary->intersected_domain[cell->index()] = 2;
> + }
> +}
> +
> +const MeshFunction<uint> & MeshIntersection::overlapped_domain() const
> +{
> + return _overlapped_domain->intersected_domain;
> +}
> +
> +const MeshFunction<uint> & MeshIntersection::overlapping_domain() const
> +{
> + return _overlapping_domain->intersected_domain;
> +}
> +
> +const MeshFunction<uint> & MeshIntersection::overlapped_boundary() const
> +{
> + return _overlapped_boundary->intersected_domain;
> +}
> +
> +const MeshFunction<uint> & MeshIntersection::overlapping_boundary() const
> +{
> + return _overlapping_boundary->intersected_domain;
> +}
>
> === added file 'nitsche_method/overlapping_meshes/MeshIntersector.h'
> --- nitsche_method/overlapping_meshes/MeshIntersector.h 1970-01-01 00:00:00 +0000
> +++ nitsche_method/overlapping_meshes/MeshIntersector.h 2010-03-16 14:39:09 +0000
> @@ -0,0 +1,227 @@
> +// =====================================================================================
> +//
> +// Copyright (C) 2010-01-15 André Massing
> +// Licensed under the GNU LGPL Version 2.1.
> +//
> +// Modified by André Massing, 2010
> +//
> +// First added: 2010-01-15
> +// Last changed: 2010-03-08
> +//
> +//Author: André Massing (am), massing@xxxxxxxxx
> +//Company: Simula Research Laboratory, Fornebu, Norway
> +//
> +//Description: =====================================================================================
> +
> +#ifndef __MESHINTERSECTOR_H
> +#define __MESHINTERSECTOR_H
> +
> +#include <vector>
> +#include <map>
> +#include <utility>
> +
> +#include <boost/shared_ptr.hpp>
> +#include <boost/shared_array.hpp>
> +
> +#include <CGAL/Exact_predicates_exact_constructions_kernel.h>
> +
> +#include <CGAL/Nef_polyhedron_3.h>
> +#include <CGAL/Polyhedron_3.h>
> +#include <CGAL/Gmpz.h>
> +#include <CGAL/Homogeneous.h>
> +
> +#include <CGAL/IO/Polyhedron_iostream.h>
> +#include <CGAL/IO/Nef_polyhedron_iostream_3.h>
> +
> +#include <dolfin/common/Array.h>
> +#include <dolfin/common/types.h>
> +#include <dolfin/mesh/Mesh.h>
> +#include <dolfin/mesh/MeshFunction.h>
> +#include <dolfin/mesh/Point.h>
> +
> +#include "PrimitiveTraits.h"
> +
> +namespace dolfin
> +{
> + class Mesh;
> +
> + typedef std::vector<uint> EntityList;
> + typedef std::vector<uint>::const_iterator EntityListIter;
> + typedef std::map<uint, EntityList> EntityEntitiesMap;
> + typedef std::map<uint, EntityList>::const_iterator EntityEntitiesMapIter;
> +// typedef std::map<const Mesh *,EntityEntitiesMap> MeshCutEntitiesMap;
> +
> +
> +// typedef CGAL::Homogeneous<CGAL::Gmpz> Kernel;
> + typedef CGAL::Exact_predicates_exact_constructions_kernel Kernel;
> + typedef CGAL::Polyhedron_3<Kernel> Polyhedron_3;
> + typedef CGAL::Nef_polyhedron_3<Kernel> Nef_polyhedron_3;
> + typedef Nef_polyhedron_3 Polyhedron;
> +
> + ///Defines a quadrature rule by a pair of points and weights.
> + typedef std::pair<std::vector<Point>,std::vector<double> > QuadratureRule;
> +
> + class MeshIntersection;
> +
> +
> + ///@todo:
> + ///@todo Need a function which indicates the end of the iteration.
> + ///Think about a better design after the first version is running.
> + template <int dim>
> + class MeshIntersectionIterator
> + {
> + typedef Primitive_Converter<3,CGAL::Exact_predicates_exact_constructions_kernel> PT;
> + typedef PT::Polyhedron Polyhedron;
> +
> + public:
> + MeshIntersectionIterator(const MeshIntersection & overlapping_mesh);
> +
> + /// Step to the next cut cell and computes the relevant quantities aferwards.
> + MeshIntersectionIterator<dim> & operator++();
> +
> + ///Check if we haved reached the end.
> + bool end() { return _map_pos == _map_end_pos; }
> +
> + ///Return a quadrature rule, meaning a list of points and weights.
> +// Quadrature_Rule quadrature_rule();
> +
> + ///Returns a presentation of the cut cell.
> + Polyhedron polyhedron() { return _polyhedron ; }
> +
> + ///Returns a polyhedron as the presentation of the union of the overlapping cells.
> + Polyhedron overlapping_polyhedron() { return _overlapping_polyhedron; }
> +
> + ///Returns a polyhedron as the presentation of the intersected part.
> + Polyhedron overlapped_polyhedron();
> +
> + private:
> +
> + void compute_polyhedrons();
> +
> + const MeshIntersection & _overlapping_mesh;
> +
> + EntityEntitiesMapIter _map_pos;
> + EntityEntitiesMapIter _map_end_pos;
> +
> + //Intersected MeshEntity
> + MeshEntityIterator _entity_iter;
> +
> + //MeshIterator to iterate over the intersecting cells of the other mesh.
> + MeshEntityIterator _overlapping_entity_iter;
> +
> + //Polydreon, which decodes the intersected part of the cell.
> + Polyhedron _polyhedron;
> +
> + //Polydreon, which decodes the the polyhedron, which is built up by the
> + //overlapping cells belonging to the other mesh.
> + Polyhedron _overlapping_polyhedron;
> +
> + };
> +
> + ///This class present a collection of overlapping meshes and provides
> + ///functionality to compute cell-cell, cell-facet overlaps.
> + ///@todo Improve design, it is very unflexible concerning extension to more than 2 meshes.
> + ///e.g. OverlapData should contain a map, which maps the entity_entities_map to the corresponding
> + ///mesh, which the entitylist is refering to.
> + class MeshIntersection {
> +
> + ///Helper class to store mesh and corresponding data like intersection maps
> + ///and meshfunctions.
> + struct OverlapData {
> + OverlapData(boost::shared_ptr<const Mesh> mesh);
> + boost::shared_ptr<const Mesh> mesh;
> + MeshFunction<uint> intersected_domain;
> +
> + EntityEntitiesMap entity_entities_map;
> +
> +// MeshCutEntitiesMap entity_entities_map;
> + };
> +
> + public:
> +
> + ///Constructor takes a list/vector of meshes. The order of meshes defines
> + ///also the "overlapp" priority, i.e. if 2 meshes overlap, the one who
> + //appears later in the list actually covers the later one.
> + MeshIntersection(const Mesh & mesh_1, const Mesh & mesh_2);
> +// MeshIntersection(boost::shared_ptr<const Mesh> mesh_1, boost::shared_ptr<const Mesh> mesh_2);
> +
> + ///Computes the overlap mapping. Mesh2 overlaps mesh1. Computes (for
> + ///efficient reasons) in addition the boundary overlaps and the artificial
> + ///interface.
> + void compute_overlap_map();
> +
> + //Return meshfunctions. Only first test, think about better design, for
> + //example class like OverlapData or suchlike.
> + const MeshFunction<uint> & overlapped_domain() const;
> + const MeshFunction<uint> & overlapping_domain() const;
> + const MeshFunction<uint> & overlapped_boundary() const;
> + const MeshFunction<uint> & overlapping_boundary() const;
> +
> + private:
> + template<int dim>
> + friend class MeshIntersectionIterator;
> +
> + const boost::shared_ptr<OverlapData> _overlapped_domain;
> + const boost::shared_ptr<OverlapData> _overlapping_domain;
> +
> + const boost::shared_ptr<OverlapData> _overlapped_boundary;
> + const boost::shared_ptr<OverlapData> _overlapping_boundary;
> +
> + };
> +
> +template <int dim>
> +MeshIntersectionIterator<dim>::MeshIntersectionIterator(const MeshIntersection& overlapping_mesh)
> + : _overlapping_mesh(overlapping_mesh)
> +
> +{
> + _map_pos = overlapping_mesh._overlapped_domain->entity_entities_map.begin();
> + _map_end_pos = overlapping_mesh._overlapped_domain->entity_entities_map.end();
> +
> + ///Oh herregud, please simplify datastructures!!!!
> + _entity_iter = MeshEntityIterator(*(overlapping_mesh._overlapped_domain->mesh),
> + overlapping_mesh._overlapped_domain->mesh->topology().dim());
> + _entity_iter[_map_pos->first];
> +
> + _overlapping_entity_iter = MeshEntityIterator(*(overlapping_mesh._overlapping_domain->mesh),
> + overlapping_mesh._overlapping_domain->mesh->topology().dim());
> + _overlapping_entity_iter[_map_pos->second.front()];
> +
> +}
> +
> +template <int dim>
> +MeshIntersectionIterator<dim> & MeshIntersectionIterator<dim>::operator++()
> +{
> + //Update the entity entities mapping.
> + ++_map_pos;
> + //Bit clumsy, abusing operator[] to set
> + _entity_iter[_map_pos->first];
> + _overlapping_entity_iter[_map_pos->second.front()];
> +
> + ///@todo Change this. Use lazy initialization.
> + compute_polyhedrons();
> + return *this;
> +}
> +
> +template <int dim>
> +void MeshIntersectionIterator<dim>::compute_polyhedrons()
> +{
> + Polyhedron cell_polyhedron(PT::datum(_entity_iter[_map_pos->first]));
> + _overlapping_polyhedron.clear();
> + for (EntityListIter index = _map_pos->second.begin(); index != _map_pos->second.end() ; ++index)
> + {
> + _overlapping_polyhedron += PT::datum(_overlapping_entity_iter[*index]);
> + }
> + cell_polyhedron -= _overlapping_polyhedron;
> +}
> +
> +template <int dim>
> +typename MeshIntersectionIterator<dim>::Polyhedron MeshIntersectionIterator<dim>::overlapped_polyhedron()
> +{
> + Polyhedron cell_polyhedron(PT::datum(_entity_iter[_map_pos->first]));
> + return cell_polyhedron * _overlapping_polyhedron;
> +}
> +
> +
> +} //end namespace dolfin
> +
> +#endif // ----- #ifndef __MESHINTERSECTOR_H -----
>
> === added file 'nitsche_method/overlapping_meshes/PrimitiveTraits.h'
> --- nitsche_method/overlapping_meshes/PrimitiveTraits.h 1970-01-01 00:00:00 +0000
> +++ nitsche_method/overlapping_meshes/PrimitiveTraits.h 2010-02-10 03:14:09 +0000
> @@ -0,0 +1,69 @@
> +// =====================================================================================
> +//
> +// Copyright (C) 2010-02-02 André Massing
> +// Licensed under the GNU LGPL Version 2.1.
> +//
> +// Modified by André Massing, 2010
> +//
> +// First added: 2010-02-02
> +// Last changed: 2010-02-03
> +//
> +//Author: André Massing (am), massing@xxxxxxxxx
> +//Company: Simula Research Laboratory, Fornebu, Norway
> +//
> +// =====================================================================================
> +
> +#include <dolfin/mesh/Vertex.h>
> +#include <dolfin/mesh/Point.h>
> +
> +#include <CGAL/Nef_polyhedron_3.h>
> +#include <CGAL/Point_3.h>
> +#include <CGAL/Nef_polyhedron_2.h>
> +#include <CGAL/Polyhedron_3.h>
> +
> +namespace dolfin
> +{
> +
> + template <int dim, typename Kernel > struct Primitive_Converter ;
> +
> + template<typename Kernel> struct Primitive_Converter<2,Kernel>
> + {
> + typedef CGAL::Polyhedron_3<Kernel> Polyhedron;
> + typedef CGAL::Nef_polyhedron_2<Kernel> Nef_polyhedron;
> + typedef Nef_polyhedron Datum;
> + typedef typename Kernel::Point_3 Point_3;
> +
> + static const int dim = 2;
> +
> +// static Datum datum(const MeshEntity & entity) {
> +
> +// return Datum(point);
> +// }
> + };
> +
> + template<typename Kernel> struct Primitive_Converter<3,Kernel>
> + {
> + typedef CGAL::Polyhedron_3<Kernel> Polyhedron_3;
> + typedef CGAL::Nef_polyhedron_3<Kernel> Polyhedron;
> + typedef Polyhedron Datum;
> + typedef typename Kernel::Point_3 Point_3;
> +
> + static const int dim = 3;
> + static Datum datum(const MeshEntity & entity) {
> + VertexIterator v(entity);
> + Point_3 p1(v->point());
> + ++v;
> + Point_3 p2(v->point());
> + ++v;
> + Point_3 p3(v->point());
> + ++v;
> + Point_3 p4(v->point());
> +
> + Polyhedron_3 Q;
> + Q.make_tetrahedron(p1,p2,p3,p4);
> +
> + return Datum(Q);
> + }
> + };
> +
> +} //end namespace dolfin.
>
> === added file 'nitsche_method/overlapping_meshes/SConstruct'
> --- nitsche_method/overlapping_meshes/SConstruct 1970-01-01 00:00:00 +0000
> +++ nitsche_method/overlapping_meshes/SConstruct 2010-03-16 14:39:09 +0000
> @@ -0,0 +1,13 @@
> +import os, commands
> +
> +# Get compiler from pkg-config
> +compiler = commands.getoutput('pkg-config --variable=compiler dolfin')
> +
> +# Create a SCons Environment based on the main os environment
> +env = Environment(ENV=os.environ, CXX=compiler)
> +
> +# Get compiler flags from pkg-config
> +env.ParseConfig('pkg-config --cflags --libs dolfin')
> +
> +#Program
> +env.Program('overlapping_meshes',['MeshIntersector.cpp','main.cpp'])
>
> === added file 'nitsche_method/overlapping_meshes/demo.py'
> --- nitsche_method/overlapping_meshes/demo.py 1970-01-01 00:00:00 +0000
> +++ nitsche_method/overlapping_meshes/demo.py 2010-02-10 03:14:09 +0000
> @@ -0,0 +1,35 @@
> +
> +__author__ = "Andre Massing (massing@xxxxxxxxx)"
> +__date__ = "2010-01-27 "
> +__copyright__ = "Copyright (C) 2010 Andre Massing"
> +__license__ = "GNU LGPL Version 2.1"
> +
> +from dolfin import *
> +from numpy import *
> +
> +if not has_cgal():
> + print "DOLFIN must be compiled with CGAL to run this demo."
> + exit(0)
> +
> +# Create meshes (omega0 overlapped by omega1)
> +#mesh_1 = UnitCircle(20)
> +mesh_1 = UnitSquare(10, 10)
> +
> +mesh_2 = UnitSquare(10, 10)
> +
> +# Access mesh geometry
> +x = mesh_2.coordinates()
> +
> +# Move and scale second mesh
> +#x *= 0.5
> +x += 0.5
> +
> +overlap = OverlappingMeshes(mesh_1,mesh_2)
> +overlap.compute_overlap_map()
> +
> +overlapped_domain = overlap.overlapped_domain();
> +print overlapped_domain
> +
> +p = plot(overlapped_domain, rescale=False)
> +
> +interactive()
>
> === added file 'nitsche_method/overlapping_meshes/main.cpp'
> --- nitsche_method/overlapping_meshes/main.cpp 1970-01-01 00:00:00 +0000
> +++ nitsche_method/overlapping_meshes/main.cpp 2010-03-05 23:36:16 +0000
> @@ -0,0 +1,108 @@
> +// =====================================================================================
> +//
> +// Copyright (C) 2010-01-28 André Massing
> +// Licensed under the GNU LGPL Version 2.1.
> +//
> +// Modified by André Massing, 2010
> +//
> +// First added: 2010-01-28
> +// Last changed: 2010-02-15
> +//
> +//Author: André Massing (am), massing@xxxxxxxxx
> +//Company: Simula Research Laboratory, Fornebu, Norway
> +//
> +// =====================================================================================
> +
> +#include <dolfin/mesh/dolfin_mesh.h>
> +#include <dolfin/common/types.h>
> +#include <dolfin/plot/dolfin_plot.h>
> +
> +//#include <CGAL/IO/Qt_widget_Nef_3.h>
> +//#include <qapplication.h>
> +
> +//#include <dolfin.h>
> +
> +#include "OverlappingMeshes.h"
> +
> +using namespace dolfin;
> +using dolfin::uint;
> +
> +int main ()
> +{
> +// UnitSquare mesh_1(10,10);
> +// UnitCube mesh_1(10,10,10);
> + UnitCube mesh_1(10,10,10);
> +
> +// UnitSquare mesh_2(10,10,10);
> +// UnitCube mesh_2(10,10,10);
> + UnitSphere mesh_2(10);
> +
> + MeshGeometry& geometry = mesh_2.geometry();
> +// double w = 1;
> + for (VertexIterator v(mesh_2); !v.end(); ++v)
> + {
> + double* x = geometry.x(v->index());
> +// x[0] = cos(w)*(x[0]-0.5) - sin(w)*(x[1]-0.5);
> +// x[1] = sin(w)*(x[0]-0.5) + cos(w)*(x[1]-0.5);
> + x[0] *= 0.5;
> + x[1] *= 0.5;
> + x[2] *= 0.5;
> +
> + x[0] += 0.4;
> + x[1] += 0.4;
> + x[2] += 0.4;
> + }
> +
> + MeshIntersection overlap(mesh_1,mesh_2);
> + overlap.compute_overlap_map();
> +
> + MeshIntersectionIterator<3> cell(overlap);
> +// Nef_polyhedron_3 cut_polyhedron =
> + cell.polyhedron();
> + cell.overlapped_polyhedron();
> + cell.overlapping_polyhedron();
> +
> +// Nef_polyhedron_3 overlapped_polyhedron = cell.overlapped_polyhedron();
> +// Nef_polyhedron_3 overlapping_polyhedron = cell.overlapping_polyhedron();
> +
> +
> +// uint_set cells;
> +// mesh_1.all_intersected_entities(BoundaryMesh(mesh_2),cells);
> +
> +// MeshFunction<uint> intersection(mesh_1, mesh_1.topology().dim());
> +// intersection = 0;
> +
> +// for (uint_set::const_iterator i = cells.begin(); i != cells.end(); i++)
> +// intersection[*i] = 1;
> +// plot(intersection);
> +
> +// BoundaryMesh boundary_1(mesh_1);
> +// BoundaryMesh boundary_2(mesh_2);
> +
> +// cells.clear();
> +// boundary_1.all_intersected_entities(boundary_2,cells);
> +
> +// MeshFunction<uint> intersection(boundary_1, boundary_1.topology().dim());
> +// intersection = 0;
> +
> +// for (uint_set::const_iterator i = cells.begin(); i != cells.end(); i++)
> +// intersection[*i] = 1;
> +// plot(intersection);
> +
> + std::cout <<"overlapped_domain:" << std::endl << overlap.overlapped_domain().str(true);
> + plot(overlap.overlapped_domain());
> +
> + std::cout <<"overlapped_boundary:" << std::endl << overlap.overlapped_boundary().str(true);
> + plot(overlap.overlapped_boundary());
> +
> + std::cout <<"overlapping_boundary :" << std::endl << overlap.overlapping_boundary().str(true);
> + plot(overlap.overlapping_boundary());
> +
> + return 0;
> +// QApplication a(argc, argv);
> +// CGAL::Qt_widget_Nef_3<Nef_polyhedron_3>* w = new CGAL::Qt_widget_Nef_3<Nef_polyhedron_3>(cut_polyhedron);
> +// a.setMainWidget(w);
> +// w->show();
> +// return a.exec();
> +
> +}
>
> === added directory 'nitsche_method/polyhedrons'
> === added file 'nitsche_method/polyhedrons/SConstruct'
> --- nitsche_method/polyhedrons/SConstruct 1970-01-01 00:00:00 +0000
> +++ nitsche_method/polyhedrons/SConstruct 2010-03-25 23:33:08 +0000
> @@ -0,0 +1,17 @@
> +import os, commands
> +
> +# Get compiler from pkg-config
> +compiler = commands.getoutput('pkg-config --variable=compiler dolfin')
> +
> +# Create a SCons Environment based on the main os environment
> +env = Environment(ENV=os.environ, CXX=compiler)
> +env.Append(CPPPATH = '/usr/include/qt4/')
> +env.Append(CPPPATH = '/usr/include/qt4/Qt')
> +env.Append(CPPPATH = '/usr/include/qt4/QtGui')
> +
> +env.Append(CXXFLAGS = '-ggdb')
> +
> +# Get compiler flags from pkg-config
> +env.ParseConfig('pkg-config --cflags --libs dolfin')
> +
> +env.Program(['intersect_polyhedrons.cpp'])
>
> === added file 'nitsche_method/polyhedrons/intersect_polyhedrons.cpp'
> --- nitsche_method/polyhedrons/intersect_polyhedrons.cpp 1970-01-01 00:00:00 +0000
> +++ nitsche_method/polyhedrons/intersect_polyhedrons.cpp 2010-02-10 03:14:09 +0000
> @@ -0,0 +1,118 @@
> +#include <CGAL/Simple_cartesian.h>
> +#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
> +#include <CGAL/Exact_predicates_exact_constructions_kernel.h>
> +
> +#include <CGAL/Nef_polyhedron_3.h>
> +#include <CGAL/Nef_polyhedron_2.h>
> +#include <CGAL/Polyhedron_3.h>
> +#include <CGAL/Gmpz.h>
> +#include <CGAL/Homogeneous.h>
> +
> +#include <CGAL/IO/Polyhedron_iostream.h>
> +#include <CGAL/IO/Nef_polyhedron_iostream_3.h>
> +
> +//#include <CGAL/IO/Qt_widget_Nef_3.h>
> +//#include <qapplication.h>
> +
> +#include <iostream>
> +
> +#include <vector>
> +#include <list>
> +
> +//typedef CGAL::Simple_cartesian<double> K;
> +//typedef CGAL::Simple_cartesian<double> K;
> +//typedef CGAL::Exact_predicates_inexact_constructions_kernel<CGAL::Gmpz> K;
> +
> +//typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
> +typedef CGAL::Exact_predicates_exact_constructions_kernel K;
> +typedef CGAL::Polyhedron_3<K> Polyhedron;
> +typedef CGAL::Nef_polyhedron_3<K> Nef_polyhedron;
> +
> +//typedef CGAL::Homogeneous<CGAL::Gmpz> K;
> +//typedef CGAL::Nef_polyhedron_3<K> Nef_polyhedron_3;
> +
> +typedef K::Tetrahedron_3 Tetrahedron;
> +typedef K::Triangle_3 Triangle;
> +typedef K::Segment_3 Segment;
> +typedef K::Point_3 Point;
> +
> +typedef Polyhedron::Vertex_iterator Vertex_iterator;
> +typedef std::list<Triangle>::iterator Tri_iterator;
> +
> +int main()
> +{
> +
> + Point a(1.0, 0.0, 0.0);
> + Point b(0.0, 1.0, 0.0);
> + Point c(0.0, 0.0, 1.0);
> + Point d(0.0, 0.0, 0.0);
> +
> + Point A(0.1, 0.0, 0.0);
> + Point B(0.0, 0.1, 0.0);
> + Point C(0.0, 0.0, 2.0);
> + Point D(0.0, 0.0, 0.0);
> +
> + //Build a point tetrahedron
> + Polyhedron P;
> + P.make_tetrahedron();
> +
> + Polyhedron Q;
> + Q.make_tetrahedron(a, b, c, d);
> + Nef_polyhedron NQ(Q);
> +
> + Polyhedron R;
> + R.make_tetrahedron(A, B, C, D);
> + Nef_polyhedron NR(R);
> +
> + Nef_polyhedron NQ_NR = NQ * NR;
> +
> + Polyhedron F;
> + F.make_triangle(A,B,C);
> + Nef_polyhedron NF(F);
> +
> + std::cout << "Exact_predicates_exact_constructions_kernel + SNC_indexed_items"
> + << std::endl
> + << " allows efficient handling of input "
> + "using floating point coordinates"
> + << std::endl;
> +
> + Nef_polyhedron N = NF * NQ;
> +// Nef_polyhedron N = NF;
> + std ::cout <<"Nef_Polyhedron Tetraeder: " <<std::endl;
> + std::cout << NQ ;
> + std ::cout <<"Nef_Polyhedron Triangle: " <<std::endl;
> + std::cout << NF ;
> + std ::cout <<"Nef_Polyhedron Intersection: " <<std::endl;
> + std::cout << N ;
> +
> +
> +// if(N.is_simple()) {
> +// N.convert_to_polyhedron(P);
> +// std::cout << P;
> +// }
> +// else {
> +// std::cout << N;
> +// }
> +
> +
> +// CGAL::set_ascii_mode( std::cout);
> +// for (Vertex_iterator v = P.vertices_begin(); v != P.vertices_end(); ++v)
> +// std::cout << v->point() << std::endl;
> +// for (Vertex_iterator v = Q.vertices_begin(); v != Q.vertices_end(); ++v)
> +// std::cout << v->point() << std::endl;
> +
> +// if (Q.is_closed())
> +
> +// Nef_polyhedron_3 NF(F.vertices_begin(),F.vertices_end());
> +
> +// Nef_polyhedron_3 N = NF * NQ;
> +// std::cout << N;
> +
> + return 0;
> +
> +// QApplication a(argc, argv);
> +// CGAL::Qt_widget_Nef_3<Nef_polyhedron_3>* w = new CGAL::Qt_widget_Nef_3<Nef_polyhedron_3>(NQ_NR);
> +// a.setMainWidget(w);
> +// w->show();
> +// return a.exec();
> +}
>
> # Begin bundle
> IyBCYXphYXIgcmV2aXNpb24gYnVuZGxlIHY0CiMKQlpoOTFBWSZTWUuYf2UAqZL/lH/2YDJ/////
> /////v////8AIAgAYLbePKHuz1mwfd67Wt7tOF9vuyvp8g6pfe5fO5xd727EVSXQfT3b1zts7bN2
> 599TPbdQ8++bzkPd7qtsU++7u2mvl67stt2HkHtBbffcON7mRRC+hdve5thsACjXvY6Nj5effc23
> z29967rcBfb3nuB66Hnn3sGPDF3du7W7TtYbn14Hb77aa+3zxlh7fcvj7ODbq+Iz733xrZ90u+8s
> G1dXTzZT1682PHvOdam93mvHXa2X21Wu5vdm7nG2Xbrdw7rEpL6YvsYo+mRfb3DzjWsk9vcYFUnQ
> myoR7e47sSVae+I5jqUusDlbrWSHzIAUEHuegda93u9a13YzvY6Oh8Qc+82667vl9n195O973Q++
> 2q7bVma1dLifD7456r53NBKJV1rB99u+xqfS31LSueuphtV1x3Hpe7jozrzr3L6ds5vffeTqfPcV
> e8Zeu5gteTWWq73d6ttmprfeH3ffbtWpdZRZmMzFStarTbsN1bRkCzNL7jna92JXGunCKF0tlCts
> rbXXOpYZBbLs7Q3t6A9eV2oYW2bWqqO2u1vfdT3aN732+9jVpRteH1t7XmaG2rVlrbWtk0VQC2H2
> 7h924jPnnhJECAAII0AEAmk9CnlMymUYTyniIAGhkGRo0AJQQwoCIgE0FR400o9T2kJ6mmyRoep5
> QAAAGgAABKYUSaKT0mjNJkaBojTTJNNGj1GJozUD1ANGgGgDJ6gAJNJIQhGTJpoCJsqfpMU9Mmqn
> 6HqKbT9U0RpvVPU9T1NNGnlHlAekD0gBEkQFGg0maJiaNMpiGjIyBM0mEammnpU/J6mSm2qep6aj
> GSaGAiSIIyaEyJppNoTEnqbBEKeU/1U/Snij1Bp6h6j1NqNA0AAHqDTaI9PgoqzssqJGCkiDUFBk
> RJFESiCKEgqodgwSooD42AKLeACdoEQAaAZdfSRAWrf/tVb//Lq1mEJCsEREZBYLFRWAlt/o/u/0
> aJqMFBZHk1qcODrrLpA1JCMUSISMRkgqQiJIsiMCAkHoAJD76fej972W7VDVtyv9v+0wnKJH8RsV
> f4o1IskUwSwREAUAYqoSIiyAiMYiDEUBViiqMYI/20qxT9v5d+oaO9K9e/hbc0Igwn7Z2pJZmGH5
> W0LNx5378pP+SXO/2ew7f6/k8650vzaN5j/f+jqIlteUF01xvnL+ysjBf1v/ccZSGnqqotY4rP7N
> jSK64PpXCwTCxVKqp3RKmbHQ1r6Nj07bGxqCfTeWtQTCxVlReuORnfqhtlJ66wPKNvzsSfxd6CG9
> v6QQxP0cBXuz9TgdiuWrhgp58D9HWikNn/z+js20w6hrbU481quJTq4ydRJcGf+3vxrmlAL5ijn2
> 2WSN+zyrilcwUKcc4eH9//eFTLW6uwjIhe8TaLjvV/bSqvHSjt5Y3rrv6mmg864CDqJA3+efoMJE
> cv7MUV9ub/b/rc0Qrab2q88qgb9+3s/08T589FsKnxtP7Wu7Pl3dnBqOfZ460jGLNkOOfGc1FVPR
> S+oYfb3UmhFJUUhPQMhJvrXfqT96YqqqiH30u297wcQrWHuYfx70DQ9HOlIaFj8EmkFmbNhsMgdG
> Yq5tg2mc7POztEOXFqHHFy/jmzNkvhVFECqqgrqJwCOb2Z/G+iXcoj3fz3bA+/Ky+Vdc2V82DCg6
> JtewTGqW04KLEfn0Ztut40Ool9t9+uurH1Tjb2nidk40yIXH0bCoPL6FoO3W+LENippPlKEYvC2E
> WpSqaxwErVh01k3LLnQonLxQ5M9Ts+CQ9ZteHnlz4e3s+51r9reH33NtLnJNLwN2Ot0aZ1XdM9Te
> 5vfg6Pa+DFrJx716CURZz/ZnGrou1umDnTpon7E705IvG9ns8vd6/e/DZ7ex6d8CkbLLmMmH2wR7
> 6+49uEH8LsM3hB/oQyzjVPSy4p0NuTU1RS5wYfQOk/xPgX8V/xL2Pm8z1XOzbrnZ8bhR77pUJZU4
> p3/D0N60uqYurqcNa6eaw8PrZh8EajK7Zaf3fp+J8VQXp9ahq3t+eGV4r49rEPdgvZj0UKJkReLK
> rCBNpNM82HaJt5vKEQcZB7/5tcXAWvUhj66h6g2ZHv3/Gw9PmuTPSX0yrSSGtffdy6b5MmZOMDPZ
> f4RhIs1dRv641SMbRGhxM8A4sxAkFYVXtA0ZK4aCZRa7HSHlSUyIgr5U+3cpnlZJPUnLViIbdeGB
> kOu4lsiyOjDBA/V1U0cNADhN5FAD0+FUUJ2DqehHgSnCaE6rYI4nR3V0c7DNPTJ5DrLU35XTO7Rx
> v5X2bw6melONrBbDoWV7OqzePGo2trkojIdkVsqoOKtglJ6MGjTj/ogIzbbgphLCE5tCwiXAgqcI
> kA44tY0AxOu1d364zw1fO7Pb1UwYKaNtlOegxOp5nbEbWu7fhuveHKEwoxZvCJPNPCTsuGZ5eOoY
> xaDybGvq6/su9I9aSQBWZbfglALl7XSwg4qmCvby25GoYIdrDcy+kftYb92/Iy1x1DhKyaQrCw+l
> zfvz3e44nLfdIph1esq3ZzOTP9m3d7bpabmUKx5PcRJDELnE8X/tLZoQlSMEYQghTFW1RTxOuTw7
> EH7w64zzdnc9IfT9O2OWVTdWv4fX+/rSJVNQpXasjNw5lRsfPvf24iP0Vh/Xcf0vzLfTvXth52tt
> CrQoHxENKzOuNFkIxgiUYCiiqGzAxmmUSKIk+T6tzNZmB9ds2+UdIVPoXM9zcMgexDgO/KM6bJEv
> MoXWrnKXrPxl4yY+NQvocvGam8Tss266ZlO22ewYfwJt42dOfzzr30r167bcEtX1Ocrod/8M2Djs
> KLukNIVO8tMiT+F9DJrTSMB6d1MXXxvmHoG0l2eXfo9LDQL52XyaybvWwdtYNV6dVNumtCE1pndC
> d7Ces9NDpz5Xg7MwX7GbbU73Z6m9222pqOVFJO5q7HPBXGxAmx8i4dmUGyZDRc6cOOcxKQMcXHGJ
> DZvwSgjBU6NS+Vw0ynV+Xd5/z5d0nH3ZgzAnQV6GMA8wVwKYKsFwM9RvH3PCyNNuO4fYnRzsulTh
> 1E07uzoQfuvqPlf1yb8snH000zj7hu6cqvGrqWlFY03zIpiUuqQRkwt2gxFluC01p016t7ScioBQ
> hSmpWmGGDKKdYOYOzv3/h/3u9htBr1A6mhtavRCMjxsv3M6qL6r4OPLjzK+qsq5wekQ2nXHy9lah
> dZjGobbz9vU+lOp6MOSTkhuckoEo7NRjDGQrAFYICKSe7bbGr1O9OXLnxv2cdetcHPHDcY0ad7Wd
> aq3LNBzYgGV+lU7BFVITrvBhYJIEPrYiOBHjFCprAOMU4wROUXAgbQcYDyg4QHKKbRMo5RArGlHC
> JzhcKUR2h+kP1kZPCHgdDHMNbmn2FYT3cEJ3AwEYooiiqsEEkWa/D1dpt12vlep67ee9X0bbPhq0
> c5I7b5vsFnId048N+HD5ZciBopVJrRAuCpU3ERqioU00Xq7LG1I5VZfKjqZV1WqbIpWoxpsmWwPF
> zmQrUyVSLpXzMO65w9PnWcRSkTVYGaszptLDEM/HGnbY22wuUEePVzwhnQxIerSknawU/DKddTaC
> eZvCJhGCScYjUCVKDdKNvD5XURRm7sxJ6MeDeWW2St0dEmRzYIdsmgog3xrxpRQ18FrDUKjGFVJV
> l/N+4pUVRt706fMqp7fljb8mxkt5L/hZLQk5PZzsmsS7HDtopUFRFdvVlfC6C/rhqTBga/cPTVNC
> +zqoGvmxcY3IJDhgd+yrpuiFWo4Tgq5z11A65cu3BAAF1z0+WgbVVXKkpvEGppGpcVVRjAMHHW6s
> kCHw6ZllRMYhCAqjBVUVSRZWNLCDGosrICigsUUkJUAVZRRRISFSUUUQgWDRCKoAVD61ZKhIkBYR
> B64iFog2gCyKDdqhlFFYLEGVBEGIFEgKLAFCMGQ7Q5YSix2ysFRINcdCh73b2Zd/frMO+2HA9f17
> Xlz5c63W++/PjXHTsiEIB2CqqqhEQiIAKBIqixBWEoohNAkQLMZA5emktAbQVAqDh8mFivq4/i5A
> OsM/DNPop33ocAN8XYODUSEg85KwHA0BE+msb7DbFPuSvAful6nBWJAKwhpQQAqBPax6MHfD/PSP
> /f2KsSr0O7w7e7v8Ed3hWCvs7T29upk8pnnNAWCbzdqhBqFCfLXBkFDc0YGBBYZqDdoIGISVqsUV
> SSBuwN5sWBKRNx31Q5BEMh+ScLhuIjMPR7VXIH/SM/Wzc6yeW9JoJpBZDrQMsGPlmk+W2znWwJsn
> 1DdTbuKGt2BoQRnaYNYkb01GKZIGQWqOX1IcnfhsaISKvBwWLCYOSbnmwDBhpCig1Bd8JsYFeTNk
> kUi5xSSadDNSKtQZBCRGU0nuPqPX/6Pj+s536bECgusxzMnHwF7Y8z6jt/8w+9AdYOwDMMCsMCoH
> hngfrM5IcpDpBnIBD2llLxTQmIBoWA/+/22tgbvvd8DpAPWb4vBAkUJAdxe1+vywRh6T/f+YlIYT
> bfc70WKYNjB7ev59Dbc5PUhwhSZ7Ify5OL9Y0OobgL/Pchxg8Uc/KprF1LLlZyufgDJYtYzlXouH
> 7IZ9ldcd4V3Vjlltk1axYVF2A7DDE0h1A8hfdsG69oXzDaSf3n+bwknAenefhzVjsoXFMLEyCrJa
> ime7H+W94E1kknE7Q9KeXruw9WZeLRbTc8jHs+9E4ftH5E5s/mv5VqI4cMAHE/UfkSb7QJC5n3sz
> SKAdgokQIiiHwJEgPvFA3A8LUOmrFrJJJJsGksPEm/Q3SQNz0Obc98480syUfEPORV6f2pDgGwVv
> fTHh36CGScVMEM4hTqajP9RgqhHgtU1Qs5L7F7g8QYKlm2AVncC3oHh3mm5667YnBH7/c93APJq7
> urzcb4B2xftIiHqiO+RHRDpB6UzSUlgaSqenpBF5/+DLuUDgsgL9lKlQCRQIQCRWQgqgskFCLBRQ
> EZCDICAEhIc4KvMzopWQkWDAPagv/YiDj1OJLv9jr1eHte+SIX/lGQ5AqlJpFUvdGDyhRhoSslRq
> 0XGrK2YlISbL9AswDx7z6iv7JAItIDhubdcK3fXp4/bij+0a9CVEidpw7bbaXM+xlUSII1kXHHHk
> 5aNs56LuwtzsY3erFaTBCYjstHQtq2rJqRRdMECpnDC0MkNMVkVaTYYhQp2dWMqnHF4N7aK2HLXO
> svTrzOambQtPSNLhZU0lmuLOsd56wK9dS4rd8dVvOlc8PxHsGR5/hmKF61GEjJkXZaQMM9JQlA0k
> 8aqpRnUpZ11O2++ePgBzSDnrqiQpoVcIWwgC+FbM33u7M4PaFlllKULaEtprzb07jz9uVDwV01vj
> TKvXVWlVtLXjq6a2tzU7sFI9rq+l83S7q264DZ2IXuz6WLHFwTW3QZ7S+6cXlFnvnUxT5nTiLFBQ
> zSpUW5ecbb5hXVMqY02PDGMsY1w9eWhWgcnFwlp/SjBmXw7eZvWsWvxa7AzyMc8yFJkYzpostqy2
> t7vGl2y9mI04+XdIUJ1hwZonikl1zO/GJLMRszODKFLWmHetqioS1wqqGlkamNarRgGe6V6dancS
> WrE6dIWcNWzOqQpUs0mWxGVWNi4hc1TYWKYlcrmGEDNtpYVNKguGRjO2cLxWCElQxyxtppRaGRln
> e/QbnvH+B5086dyalgsOLQ0OXezsPHj3f2vkPGfdP1NerPX6vq+I9829ff/KJnwcgDEAoPZQeUbA
> hcGOAUampZ+cdMhtE3RicJWIwUy01aTBEgYyAYbRNB3aWaHRptljPntAJzYIULDLHlsGEbOzzc3h
> 9lrFFX60/F3viPUfsf+L9j9lbF56nJXeiHXmuyD6jxrNoxfdR6C7Ll05WK9jT+B1gh25EzmSEGB5
> cm/Z/I7BGQtDTzUSdQxVU/933Q1eIfaEjA7y/q04++AZ/VPO9urdNPyPybJkmwBhhEulODSWALTq
> ed7JO4OEjy8m3fY6niuaWdnktDc7bbbnd3HHKLDj3oZuJTmBlrLkluBxXFbxYOVmWZyQcORwVqnn
> jy6JhxioAwIQZrOQDKWhrqXgSW8SwgoQZ05HbGHsESSZIR4i8ha9BL+JI+LuhfiUXNchu2CSoULs
> D+Ac3/dyP3zwPhQCcgmClLBu81s/2M5hPpChvnXw8T534NJMIXnfl16o8EAMnYQon3h9I01QhSnq
> 1N8LMH/PYKFwoVg5YKU62dfl/PwXH1BNvn7EXCsYppqRpBqEgrh01/MSioxaofqC0JTZfqb7/8ZL
> 3R5z/jlttuiZkWZLbbb3kyVmS2219QTQBMfXQFsBzkXxMce3Wyvi/DarayRcw1wzDRRRo0c/yHmD
> kbPwPZwatOBH/2TbeIUf9kfGDyoJcuO3tsf+AG4BAb/r5EeG26wBdW5oMCzHZxpO8f3hFQtjHbgy
> MQXF/BtOnA0cM+oiJEQRIiREiIiAxIxBiREiJERBPj0j2nPasQOGZnwKO85ev09xO8G0FKRgkJRR
> QLIlFmYZDju2wmNd7wUB/xjBiDNT0DiSlxvxoQlPO7NbhtkNDPDlLD0CiiiiiwYxjGMRFFFFDmnm
> O/z7Hx/f6OrjYs3npmGTJmSSSSQmTG2x48mQmTJkknOeqaA+Tyh54EzMd3qzRzBAuBYqxpS9SVI4
> crFQbKiooAlIoSbUiVlVVIjjP6vxu9PrAOGCtpCEHNHaSnXt2HCyTrPo3hnhfDUzq2a/FqWmZ1Zr
> JHHPyvWnxuSSSfe8Xbu6BEZ7Pvyeo/V447bAR032AJJ+9CYmZKuEHp++BIrBaMnNlbXA6q14lvFu
> h7T7hWC9goo9sC7J7ux9Ga7blMb7tZmV3ZoRYvZRlxeBUM4xIgspqPXwlf4BQUL9WbOA+reQe7un
> k0JIhA0ucgPVIkQlYzejbhUYsZmiEvd8WLyfyNtpNtMC5rWaMJk93x6NYW2l1mKVg4KSoJxRmqSm
> l+R2kJK+7bVQGGBQI5jrn540aYI8Jw+aQ632zVO6+ijkrk/Kzj08BrH0+zr68NH58/vP/NFJWgtR
> RknbCij8+MGlzY7MKQ11S6EJkpfmyvhfadJLax7WeRCFvZW/cP1Byy6I17EjN2DxCSDg4YuKWbeN
> 7wiWta0LiX3FKQUxPGCebySuLLDNj4CrQPiN4yMYjBpXtDjF0kK6tqkfkipzni9+eby6Oy7s0HrB
> 1sXl6b3JrZz1xH6r0Cf1WZVQ/gXzftD207oG21cKwo4v7V1fnp5h7CdrU06hjWXVVykyuXcer8sN
> 2eGDXlbGGtRn1VO0Mi2F1SEo2dWFeMuszjRFrMxwuYla/tbpbNr3zDQJ5Bg33F3DjUqgf9NhzRrR
> 6Dt4eoDD5MKnJ4edu2iIPPfA2HqAsSwPEQPc6VfBTg4EEF+mzTWVWtk2bJ1xy9pV88t2rjfG0i9k
> oJo96g2bIM4BHTj17+T6nvoehEKTgvdVNeGcsgPwptfiDsmfZhY5ZPvIOiZzA4507bViGgVzIvkt
> wakIKthBGWT8s2nFngzwZ4M8yIbILGTc8q7C7eNZNyOlpJMirHte8EcUWwbvH8xq7Qr34MXMR534
> wnQQ1bZyjHQggyEssXSN1A8iAKLeTMNgYb33YMzcet54X+WXNXi4DF9V72OffV+bf0zUeqpfsuXA
> 63Jx9h2hk8A/3E2k1+rd2c6EnwGKLERRQgqCERWCkBSCgx2+BDp+r+o7uHi6WQf2nIsyDjPIpfB8
> UashZ7yPjDDWUJeZEHr0rwJU+WuUSqgilUxK2z3fCJhZtpd/BGu4HCPwwJLDgsai7zwh6dk1vZzs
> Q5Bi7wHSrK2lCpo+6GUth7KVJWgMNsOmMrbCF/mguKJwHiw4wehOy6pH5p7yyo+LblQPsM0uHUXO
> JrKSiosWzyiMRrWMNqMRt+Lttttttttttttttttttttttttt22222hbZC2yFtkLbLbbbbbbbW222
> 2222222222222222222222222222222222222222222222222222222222222222222222222222
> 2222222222222222222222222222222222222222222227bbbbeO/s1x1ek8eT5Wc06Ont98eYgX
> J3IqwZ1Y1WEh6j3TpIiUfZ7qG9I+hQw/i8HWMVwXQ9I5H3fM+i46QXsKrkV6PiayvznQs3PbmiZ5
> vIDwzeZuxM/Aa5W0EN2jzrZzcGlppk58AzhYK8iy8/LDePgcS52tu+Wjn4uNhaRwUHU0C1WZ+GvW
> BVg+zUzFsr2IZdUFUhVVbK12tkxOuXposySxTS51E/CtODhbTK6MCaebMmDMhnLzGsv/VMuhXdrf
> lrz9/O2rZqvhd9BdOwQlWwhLQ3AytXAvuQso7xfaPSJkx2FQqwyBjrhmQTqZaGd6hlE4Rm7WmVGm
> IfsiOietv2u3mOeCHtHq7jNYxqhYWiNHlznuO08KI9TWOyHcOd+noDbzaHDcHwNndier6Zj2Lh+r
> zNXoWky7w0Bo306m6L3EawxhzXvIwvk1xtmfs97nAqWKZXSGAY23erI0RT2/M73zn1MO3qd+XM8t
> E5m23F49PQX1fgq8lVZ0sTDKbuO7w7td8V3AhLcvhqPwchke8K0uIHliEz3GdSMohcGfPHNEZcW5
> A3AqeEMoVMzlhakSkgpVLfkYisi6PfmM7oPozaHvhLgq5Y/G/BHt7xOUnYOkAlG1j29MMukA9yp2
> 9KH6LdVVzWytS8/OFiktKgvIPENkwyYdFtHtj48YEIQgQeEL4+DP18LK4jloxsAMXEr2nC9yCTAx
> xC0D4sDDA3AQARInU4G9fIvI2XOF2NWccWJNfqy0q2D40/1RFxHS/GxtxtuTHkRzgD0BzKtBRkyz
> aebVnAopUrAD5MUBI5CvGSJzss9FVM0PgVKmWhBTphmMRp2RISExIDEIEVBXQQ2ky8Lcm1bBVNad
> qiqDHzetE6hACc4ezlnErw9r+uT19rR51kEy2uFSqd3FN5RWhfxZxorMvTZ42wcZla2vWMq7LbD6
> +QoRpaNQVnC5nLrDuZbFYcd0ZkZnFwoxs59DxH0ySj4hW1lbfEbEdEn1s/GNnVGbyX9taXPoZ34y
> 6OVqMuNg/e/K8flXZuYXr89S29ytqROxu1T42PpS1AtLDSxRa4FTr3jYyuDE+3HA1MdYz4GpYCED
> 0nJHvAQW95DL6fOMzSQhKwRHl2lGNQkLIGGaGbBq+YeMajckwcc1ZCydl2hoqJBghLh9H9NLJcRz
> wGOuV0DedAzi2cYemEjhEMLjWo9gBElBlfs10ezxadhePSJU/xqQ8tIXtF4wTQ5+b2UoweL1ylzr
> doBSZuPDpj2Ii32pEX7vA0rtflA/Fnx42tejhS1FkGEh7rZrtuO9k9oSRAmSfGJIKCxTo0ZB/FOl
> 5At+BMSCbA6YIGZyIe9x16RIEL+mROcEYH0ec/T7ffXtfgm/3pw1x9/6llZ8snt89NvNrc5vWNS+
> 2y6tlJJl9ackxDRb5NKa+wVK2vWs0ouM6058lpzDq3TfmOx5zzBJ/AgoKz/Aqp4/Kl1bcMzM0B/6
> GDZR+9D4e7rzo+07PQTzwvgiEu/xI1o9bvmGhZvg3qF4FwaZOOESnffBWqYUKSTGs5VTGmKjIQ9Z
> SNEziRw+9OuRfd+n4GhPk77paWH7PDb6eWxm90yu7iVlVmMzfbfReXvTL7+cz7jpZi2oikWYa1jf
> jnnPFFPHxGXvphS5bLiBkIKpEYGpAWAeEk7OYdn8NOw8E15Y3pjn1Objtmrd99tMW96Yqwbg/Nfx
> 0mtb30nRaWtKROYevRu9NRl3Pd8emPDe3h8/swdy57ddVmtfJci+u+MYr4mTR7bv3dOWPeY3UyvW
> 5ydakbOdanrnfeTY6cxvnan3F+Zx3OlOu+uY1qOJ13U6rWG1mqip1Tu9TqqiaaqznVNnOMOrPNC1
> Is0LqhZoWaFmlqqqqmhZoWaFmhZoqmdnFmhZoXVVTxrVPNZqoqdVU6zhFXBS4KWcqqyZWMOQuiFj
> UNCsaWTS0Ygyy4ZZ0sRpolZws0uppqF1STS1hYpZpcZmlmlzlZpafDPUxrGKyuqeaxnOMY0udaKW
> iF0QuiF1MrFUsyTTNNLUlYfCvWVqdU1PS1NLU0udLlcZrU1qZqqWaMxJNVitTWZqqxVVVZWsriWh
> tEaNRTS2MVhaisvFa1iW0VAtQlVSzNLItTVTVRVRVRVRVay+daqJmqmaiqiqiq4553qTi8W3nnjn
> m2219JNa1bvXPHPOrbbdc75ipqpqqqqqqqaqaqFxFZimaCqzjOomn0Q7408Tqq1NVqpHipWI1UzN
> Nlyay+dRE51E1nUS+SnxnOtROZV9POcROqqKqqitYp1qqRnxnUTQ3Nbb8cuaKKKKKKLFrmtt+OXP
> jm3NZ1E1Kqz4zqJilZ8Z1E1NKKoqiqKoqiqKrJvfH2NJ2z7TzNGDz7BLo9WQwGBTKGAlJL3gAsbd
> b3c1NTUnJtTU1K5Vs1ddZXXWV0roaqKqsLVUVWiiqZmrC1VFVoqiiqoiNJkwSYMQuSDRROFog0ZM
> MuDA5rNOuiiCiZWiRyijGk0jJpVNJkfCZSMKSlI5CpQqttVVRRVVbaqqvOkFILhMppH08YkahpJe
> B6GgmRpMvA8DQNWnVlWih6KKKKqoUooop1oooorRVVVVVVJVVVVVVFUrDUNVgiAicKiIggJyIhv8
> C+uGvV/r27Pj07ayzxb8uPQz5rWBATpkQTjEBhhx2wtlIMKyYWTIjAyTBEogkKURlGEwksC4BkMC
> mSRYAVAKkKhIREFCBC5LCSskBEINqtYCx+p/mIQfUEaEPVFyYAGhFSQVh/SkPrfa5UYjlCfAa6tx
> CBARkB9oJPeRh/GJ4JCnxSk28LIbsDYpSlKH+MA3QD8RkHlj7518bE4UFGpmd1Ws/yDBZo82CVPj
> hphhhceUfRBj61UikYb4FFtPmw9IbBlALwvOvB3HECEWVEWzYHx/a9Zfhtz9+XHHohFVVVVySw2Q
> 9NZXRuywEkmRVVRF6Q5xV2OjzOOOLbSlKUpeAT8c20Zm34+viW9s1opz5sPOqgkNuXS29UOSIqkk
> hYkrozWFE1TPq61kxa+zyErwuuK3okz5fvvLq+UvXf/V4Vdl8Hgcxo5EZQ55Vdjm5VNy6nG/+czn
> wcTN+JgCqGQtLETv5kOgihpUESH5fbYKAeADB9F8j6H0iILxkXe2aqKFHTfYcdEStKFCE6QWqJky
> BM/l/dE5ImOsEfxsnOT09hY0AgZHmSFiczU1MlckIARPP0zrpylKGzq2ttttpbars0wp+kPTA7jV
> SFpj0oKD2mbpCCYmIY8wbHHvHqvO0AsP7PD1ZAG0A946W0Im+Ef0MSwwUbybCAvYEmXbau5Af3pz
> yz4LvGOJmdp3n5y7du3Ln0H0H6jm3jagnKh/gdNfy6HJg28UA5CwWbznPISjm6zNGoITGwOISbQn
> xM0ZrMOLl2LpEvkmno4rSS5PNhc0XRCEPP2ymwGwgIhEZERBgiZ/r0X3EspshH8Srp6B3zAuuq9k
> sHrJvRxJimbR7OWvfyycs7mjj6jzmu/T0OZkztdzBnIsRrOAcxyJKZQhBoDsXpxxxMYEAuma1H7F
> TTnJXwQFESTgWkjqIAVFdvNiqiBiboP/QRNlIsKj5+XLj51iouu8rgNEsPm3E5mwveDB/Mf3H2HD
> x4n0tHljV3RP5Py+VwN1p8WVR+12gwY/jBSbh1D5Tqvtl5f6jS71vLPWKTyj3DbtdOZg0dJRd4Ec
> rob9DJebWFvjaOX55rvGfhBV3mjNk+pRPBrW+W3FHLlS1K9vaSdQbA5v/DvAvCqKOnaB0bgzHn6f
> Xxp4hr50cXeFjCSJBQVFUL+3LevuPz3+fOlR9fdOf3JrULEPMC1rEkmCEaqVPSuLNkZ/7QXDBjyy
> BQO+Hj3Pw4Ss7/ugV6ye0TLD3/J6Pb+HPZeqBv7vQ+Hvqj8i0oKjcTQAVDWPZ0H6S3vKAahDbLDz
> ErD1iLERSMMNo4y/G6UxLtu0oAybJP1JuB1k/UGAh5TBXXcQP+HfdcA+UHO3U8VCeAgj2/2yPh8Z
> HgyLmSNhibfU+ecUQO7iUE2koDbM4jxcISK5z98ZgT/R+M7eu3rylWiLuwzO84PZj6/KUqZfD39e
> c3awaytGD4SCr4vKrrt8IxjGNPj7bq7kTYEQAVdLa74vH5Y0R8/Zl5R2a9rLlir+xy+GELaj81UQ
> TfnzhOre7Vrd+W/zaNtsC+XLLHk8NicMs3HOFkPF+kG68L7nrz3vatjwd4cr5w4W6vx0LiwaGd34
> 2Qvjg2zAwjuLLmC1hBfH0GzAL7/j1juJr4x1oqWkkkooqtXz74m9tqrLbF26fDu1hr1ePjr4yOJQ
> 8tYxZpytd8B3c7uhuFbM0URCKoqIbgHnIHnBBCezW3p8fV5ajGMePLBcrXJebduKmHE1PPCnbYKf
> I5FRpCDNDlO74NdDFlOxqoBxXi10ZC6/X8pHNasmGOZs6hNTbqSQumi8O3M5FdmZSE122Y4YzreE
> ZRlAZ+ztoz+GxeMLFsun5ccHwnPAqVVK4G2udTxd3Cz6z1h4+5ylCuBf7Rzs0w7muHPVAo75nDQH
> nJnhOuxQzMm/JMwR9Zhni+yfT/k7+b/btn0I/Jv6VKj8E+l5n148++Z+6LAaCwqFEQIMGNylpiEg
> w/Eo5elaDt3J3Dv+eByGKdIm9s+4L/mwhD9KaqDP6UNyBqmqRIkgzximySDVQeIHqAnAaRgIwElb
> 2vLtI8ENQ/h0QsLrDRmsZV+MP7GUFMi1WGAZ5hezT7DmPaYJ+/RaxGGTs1/G6oQiAJzQkUoQ8/zF
> NB9B78kP9I3ffchO6LpESnrKCEWEJiVRchB/oDeCS4IGSicvEfJXP0HkLK3QU2+AmRfh1PWnxGhc
> xhH794A8QP+gxRqCLIIEAIAlwhrbayZ/vBvhtvD9h+yifnDFPYvtGJF/L+rxPSmC4LU9oPQ06EhI
> XHvN/8QenBwDQjGQmfWO6Mj/WAP3aNnn/f5nd3b2uzzfsz+zQ+ECHfptGdfIyn2XrEcbFNBcEWcO
> 3jgXYooi8g4s7HGZskp1xt/XILM4+X+OisWWYfKFVgKyl/vYuQwwZiHIjVZKtSJR4w1KidnF+XOr
> YxwLWYse8Z1VV53Y1BQ/E6Bx/SBcBwAHAsObDMQXRiLCIAXIwlCWAi4AQBoCIfeVKbO58LmGL959
> 5WOP/3lZtbO9iqMvyyKTNLB7777jy+/83TuV2Y42Ru+b+u4OIU2jIBIGgMhIHxe6Se+EVEuOBu/b
> 7wmpkSZOGHX7woRKFClFeCJiFgqirBRZvBRYKKAiLD6+fvJDXRfHLLz4m4Al3HH4lnRzdjCCVGEC
> QkQ+s3j/ADkQP9iLSrB6FgMAH6h+ogn5mBNQSeAUoTBJS0niHvIeAF8xQQ7uhDDCwb5+ps33QZyD
> /QlGMY/5J/H428+6X22iqLOrjkW1dQ37zQ1YUsMipCEOG7cWjJIm3aT0EPSCHpBBCiIlCIFkGIyM
> SqGiEOp9CEHq5oPaZJQpbQoeqAoQIDovGTDg7i3UMWmhsNdK5zRMmTlCSxpGxr0HXsPjH1P8cGdm
> HSk75290mujkpovVvf+3ttmrqCHMQxHRDMMmybwM0u4gWaSyIGxuGhYiOwMgIGAQIGwDARcy4K5j
> YVcQIAJBigbueLtDv7Dlm69uHUpdaoQNzY2XVlK1ROotGJ4k9Uh0EteQhIQIfvgGAbcTpJyPzLci
> w2SnrO1Ib8KIoupO7yVTRpWaCHMFjDmRzchIER0vaqq10uHuLu7bQ5cy53kq7hakpozQN41VhGnf
> hwP64SG3lH2O3zh9PmAMBOfP5hLIp2qpn86oes2PoGOKvSL0C6pq6r8AiZZk04QyxVHsTnFsHqdQ
> 5FWA6I8Be8XARDkgPJHUMVLKEQtdILklMXouw6/cH4fw8IGJLj43NRq2lOEAbnCCk4/G2EDi1TPZ
> dB2ql3UjDtdFvCk2YbpScW7JR5Ud/WhR4u+V53oEfX28ZWTDZNb+jmFf4Hw4znnhnh3rX27Ex41w
> nU+7U4PVU1NqsWgPtCDfBy+qvHZ+a4MGlHRPokHkeeeDIoqqShRVIsBQPSG5EEEEEKVMhShJYRKg
> BRERVhGEVJBCCEERKUCUolQ0VFigoAoojIo2SSCT3a2RzxHbRpFUc0ZiKtyZl1C0kIIT4fRhJpIw
> UVAQEBAQDhhFkBwq4mGPiODDEZATK5Sr7ATTxHkHmBCgVhANku35uHJjgKQPvBA5xoSBcq5zFMTI
> s0IVifenExu5EamL20W1GiPnjWsscplgy0oYFgoCh81k7E97rOHPhi8PA3libOKPgW0La1j++oIs
> mSYMs9RoIG311iGTspNE43dVjRtvaw+5h3gyYvAgSPe7uvX7HaSOa9/gjNsf2szMzQTDDvXj/P8R
> dfSlLW/fQcNzMwPQ8fMc8iR5kyZQrr+dRGZWep6ka2a0r8EkkgyOG3OG5O1OohDsbAVgIKWggSG5
> qB49S7fHerWBcANozXggnPZt3hCDWdcUvO8+rBx3Fre5iT80mSdMMwkyxJrQ6mhhNaY0fpA7VE/Y
> oEEGDFFX32CRkeSk6Z7Abv+8/VLE/3/5lxfy/5KQAmRAiMJA5JZyasmuWRFigoMRirLDabs0GwTY
> BzNgY4GfQ9aNDQJ6HktyESBKUx0YYxNrHaCCgVorJUUUUSTxPWqaNaU3ZirorCx1Ebx7hVEmNCqG
> YbMefiHgeVWQFC6In2CADBFLI8ER+CASAJ/MCJEjFYLRClP/4L2Fu3ogFDxCJhkWYFEHwfWefH57
> fsO226D+yvO2CGssLTRk+xN8ec92oGzDhjxf1Rrpl46Tp9W/LfffM33+MJbQtsAtWkltLbAtstoW
> 0LbbS228ihE69nf6fpPyn7TmL/GczzZ0ttzqSuqdc6t9maOxe3ZNoP5/v1xtCVJ0j/MmpnKtqS97
> KvoBqAdh4pv1Ha322XU/WI952xlN4+EqzpwCBYIAq7F+GKPnka5eEWfbDO4YfT0VfJPkMSUNJkbj
> /PRRfUBCoLrEeEH9G2+AeZgTZgsAWbashjNss5JCockmIHHHbzNSByZG0NKhjOXKhphNIFUzyrDd
> SjeCvKLprwgQNP1W3Wz1rnrhoZ2wvBE/Jng7Y2+vOe53fjFVxtutVecM5ZUqtJEE61pbaGWIpZqq
> q7uauajV3VXaCDatiCDRMkoIazU00Xmc7JKfrXfpBtkU/p9oCKaF/0M4mGYQK20YQC2coVkVGGFU
> NbdBl1ES5Eys3egXAMwyf2rLCWWAOYNW1UGgCyo5pyETRGAmRgm4ADMCKhAIIAhGlGM/5yZN/Sv/
> F4qqqqqqqqqqqqqJIMIa3TMVWMeERVHSpVtHiqaVWM2tVWMcjEcg4qggYuBcVQZwQaqrGaVqqCFU
> SqoM46jnISHcBThJ3pVgLgWBRwO+lBEPrgi8Ry/eQNKx6DANO13oZJveDpLUSAYwsIWFVE4AomvR
> ERPg/CRjOGbn2cUdBSXP7uSKENpK57sOK51kmG69/AshyxBOQmbcFGmgVcQdr0g6XIvIhlEYGYQQ
> z4Og1nY3no0OcVlucKWiqukhDACEtTOsLr9DEtxtYgYFtLzHO4oYksouK2uQmASuhNZGUaJKNRrC
> BEcvzT16QNNZzSJ36ZFs0LUuSCZWUwigLrYVuMzlo4isHCnoRuGZDkLNaLo0Vic2iNAKIPz3OwQQ
> ukVhuuUs4xv56elwwg3MIiI+eFBtGCQgEQTYUzIlnUGDlqMXXHWmll2ZhtrdEd+iA6KpQnATEMHn
> nU2LGPA0LBru/agGgC4wVOCL3HIyIwh1SqhT4F6gQgXUwmikwxNWYF0LAdNIQkmhkM9kn4gKHGCH
> 8QXgghk5JEiRkSEkIiiigxjGIoKqosYxFHvNpQC7HUGzH3G3EcrSGtltqKMwwSKVELYlSjAomUlG
> VRQu1usL5heIqDeVCy4fnyZrnfcLCRhYZs5EPsEXGEbxuZAzBq4ZEBK0SyM41tWG0mECi+IaGho4
> ZGbyGZvBkafkKetT2ifgh4KTAyl7Do3U4DRL2bXSBpHvBKGcmqQi9gA4YApNp8y7JJMS3nxdALAF
> JvaN41IRGGP1+bHA4FeZPCnt8Hswd4P2AZJ+aZpfsJ96YJdLXk3nBJ17tvXxeDFXltS1NecpvLZt
> j+WZ1WGfTDwqzUUNZqXMxM1sfj+R6I6MeUPPnffhEIh2da6WOeWnBm212cqE7AMwSZJRgDiQH/Dr
> Wvc/xx7vf79rWqZLdHGuiRqqU1JOpmnu7u1u7urux31FUXdXubwwgMozMMogLiN1qGruLiddQ1YJ
> mElbK61qrY22Tt8F4niakDYbmvhOhfH0G+e8n2n3cI9xTacbqqwQ9cH5AymUZtRWRyq24SuFt/Iu
> XMrs5JspBkGKQ4Tf4JGOBCE3EtIxsQhK+HsTY6oZ6ur19nZ2eNz3zSGzmIcJKIVJNrxBzgnC8LuL
> AhfR0fF8gQIEDvGAYlIpQLht3A18vX0+7H6vQWCwBjFAYKd8UiIxoIIiMk8MiqwRUFEqqlFtFqgi
> 1FEge+6J17TAYmytRtHOoPZ3MPEwuA78xURWM5tVT0BE+vRD5lM+BNdD5enOXx+W/okiqoqb8jg7
> p5hwH0QRioQ6N5UcDE0HEuKi1FDMFWIAJUAWXSDc+313uhgCj89kBsOBhQ8FY9VikPfpodsAzlrU
> Lazaw/OC58B2BfZ9AuAebp0JqnPyOxzg18O4LhedxX4+l05QZJItfQ4mozhgYxLCusWCShzBooig
> 8yBBQFmzxul7D7DonswcccR7L2sFkYVALNU0VWiB1WRjWJcxLAs2PeZTOoOBft6YB2leaXEa4rvz
> wYvWXuWKKo7cHkPcjxca+Y8s6lVRu9UIKFhY2CxcaUoF5gkG+wSk6EicIDpO8fgn0EihBTnJRIBC
> Fsu1BYTETnOSN1IjGg6I9eouveIfVvv9A+beOsfXLu+vnsfKt8W+3njMXUb521nOzyu+7Vu/fJrh
> bWtuHu9+YmvWm26u3oRPYbnv752xJ6wOeYr4zv9OKqKnhwew6hp0PZo6cjJJc9gaKjduRio4gFIE
> L1KqEwDsiTFFx0+oiqqQIJMCZKWuZBEIuPMvL+Shwq5DIyPb39U0eBKTtztTxkHzoqqm54naNDme
> Jo6BiUYsNN0kjESBTnoDp6Akm4hqbqB1yQCkoxOHB41CeAztz25nVlbPyEvaJd4lubLQwNtqXDd4
> cZtGdAbcbnLbAZwjJJJieE2SSg1PD3IJE+hUhYbRWOD74oSIFWQN+tCRBOkDLpDQt0CFtOB4U9va
> gPQlsSokjuw41CQjOWxuWjg9hvakGMy3W4c0QkJkSa1uoa80HDpTtod32PUSpbjra6RLVCxj0Zra
> 1VBmPBRQfKuSkbpAE8OOoL6wOhG7cQ6uj9fOEg2H50c0zkluc9OtQDn6ewyuccQCZiloinNQCUcj
> mXArkpFR3WVRzmZ1SlU8xFudf1civvy3DiU4LBNwKoEcWZG5gkrdjOkUJiZU9b5KAECRyOGHUTwB
> L5LCHf0lLx+XyoRmTgu3pYrGd2B3B0vrkihEYsRYyREigKAiKCgxBSIioiogKLFkRESIMFIwQRAR
> ERJA/AA/GQn8CBm8VixiSST0gdnJ5O46TuO7o251R5aDf6vEWwJr9Zk+pQr6klzkbQbYcLwZ3swd
> WwUZPzCUOYKgLODOaz31vhtyOjdJGPMzu+X60CqaJZLnIJKB8o1+qrD1vrr1CBpplEtcaZkIRgIR
> qkc/eIQUA9h0b0Xdzi+7EGYNChQP1jCMiEZHB4bZ6+ckuNk16+M3f2fs1hN+LyfFIYj1Uq7OlUqa
> Yh2gevr7IqKRYfATxnP+K3p9Gx1tz02mbyJbYkvTt6U75a7aMsU09oiHlOvHj4R6qfTJ4OPO52o7
> rkk8XAfR0OJkjmgc6btW7eTUhwhZvwvvtK6wjg94KGtZw7QQka0Q0RLW93cLdXU3dxM5l6uLq/xi
> acyRC7VRRZYj4natUt3fG+31HhP1H1FO4e+dh75fUkw+8PVgYv2pvlVMLwO5yxtJfQn153wKJDcB
> mpBKIhTvK+FwLePj3Lt7pWqouPwwEGxqX7EUJUrTo6xJNcZN4fp+2xhcGQ8BOIhCQkJCI9lOCvQF
> G+9SHy9WWBn2hXW27bdkibDFHux+ITHD4u9MATOIDkIqh4HtMz3/YYnmc+JxueMvm/V99XmYKOpY
> tZ1+YENNnK75lbNqJsrZIBhtyHBFMcd4FVj2T9/t9UQKOCjDw2RCcE8QiB6JmMMbAninuNDShZPg
> rRQHqZvpVWuFVDx4wmxrSCIOahA+4PAYNICkHQxpFGV8hvxl2aA5SjwVuyadzsjCM2luMlSG+gID
> K4L5XYhUIxMC6y3W4ciLFpaG6qrLzsygdohhBej8maNUcdFdkydUK32VNup2Gd05HDrrNMJoVntU
> TiAOrYtAlQdAKkhrxGdo6UhDogFAoYhWXLhvtIzA9JxL8ZEmioTB1V27N6nUsZ1gX4s7uyGRUU59
> KPAUjaNjsteMO7CMKqncyPyrAM4MGzeT4m5tNINGUUAqxCMyAZDbhatSnHgbh/wSSLdtmrdt7obR
> YTYvZLWOgqrphsvb7vXco9kpkI7HnVyHSnK9KqJ6aa+/a9WRMlbhx21O05LfFjZfs9PDrMwM94NC
> fvMJIPCBiZMG7/WAkH0SDuSR3siIGB03L4Y68hRJ6C6+LAyp8fj7T4F2RMo1JA/wSiJNjBCtXuDY
> Xi5gwoqOt7QgphejAmLEqFbfjUjxlWXkEqBsd7myWjhoesYENnt+PifM5Pieh6FnMm+EXjKC/JAG
> aSUwX3JGIj5fFX8fUGF58u4CHPMjybPx7vUs+dAb3+gQe+ObM6J4eXq8IomhIGz5qzAEdohZiYTM
> SA+wM3LrseELuyUJwgYSamLY55he1qIYTQDV99NyCAPcEGpju7N+/tdYnCuqiggO3Pv+gb6PpzVy
> za39/pzzVO0wXBeKSIJxMzNZuy7lbu7q7uavWo1cXV/cIxq9JtJWxaVitrk0q52WV2yySDtNVmMs
> yIGYwBWcVu/vk655polGFC2wMalh1I3vcrnVxPj1Tw5O9lD7jgW4nE3nAF5fXAT5+Pz4wFYdc6Aa
> DiSHl07+nCnVKC7O+TBPsY6TUpdvyc1SCKSchuN7S5JKJWRadcfCAW2+HUOZfWF5G8Ks9HVL1fBk
> aXxuDXhDainTprXXrpiqhniNtQugC7aVpGjSSYPJjfqnTxnp1A40DoI+M5krwIZPVgOCA8xIrwNz
> vNGXUOTgHXjrNtd1w46Zm7UvwtvzkhtWw64YYqueLcQvFWaY9nF7HLJQpb9QloSwWNYNYPhjNkV3
> wSLq5e7XHSYVmlVaL9VAnsVemVbGJoQHthnFG0WGsYdkzMSx0iwGVdhNVhMrsmIQVDIA3NIwGDNj
> MufCxMX7D7CarCAJgQ/QgcBdogcQdnPXLkdGj1TdRw4W5hxOpoWuiGkBkSQhOaG/kJfSjDcRgpkj
> FfikGtirLyFVDIjYiaGNGe0V8lGprnzYyL4hbebmWM8wbFXJsUGMUAo5UqDWzEtNt2WlXxwElnic
> CtfOivAXhfncMISd0CmGXXXfepKvcteAHiilEEe4ielQgbDn56NMDEdaZBxCHIptXxGr3N2qqOqM
> 1erkYX2Kr6xEz593lPjxp9vB1JJnXYewHMeEz3KWbHk3g8uGjxjOP2pB10QjawMDExtCvEpheRET
> MaWLQ0IEcw0NexHLeCW2/Yr8wpgPnVHtVp4Dd5vSnNDpr7R6ZZnSuvcV8phC80ljGrwkWBKiafNy
> 9mQn5iarHyFx8GTZd8x3DMxz3i58RSmFMXq/sE5zxmNi7p9LGjynLsKqiGHx8CX+fe+L2ne9i6VF
> JKhVmYSCIuZpre7uFu7qbuzL6iSrur/EIxtkhJggskTbE7CvWrl/pPrOT3nJ7E/D+V1V1+1jDwg0
> IRGIzHH4hL58l85hcXFxeLNC4JDI3+fXFNoana5cVmSR8I824ROhSIbFXAt8GWWCIgeoVkaGMoEw
> kNKkSMooQg5emg15qMX1EWfLy+6BGNyE32D2UbJYaNJkK8ZN1+gvT8CRe/cTjwbw2TAaTwQayd1s
> iG/jsdAEnmNk2MbInBv4iOzN2gTy1HlshKIinZgs+ecpuqV4SqhFAKFzrNnDSorPbEevIV9VgpZe
> 5Coi/YnbngWrc2NClIFCGojO3JVlxsFGMzYJGqIGtVEDvYggbi8O5Xg2NUpr1Bg1wekGjfxJ5Cih
> aY/WXp0xRnUA2hYLO6SG0nlDMqMy+RnUMUBMXmt2tRcuBgBUPoo5hZfpXlogxxfNqUeAizElIqoH
> mY9XRpOs/wB32vonobua8nCoeV2p90Dz6v5fPsPBouMV6mW8nY6EsGCDDFjDqOG0ME1xe+CDCthI
> kgCWcI5TzxMJX0IF5gFmsDuYLILFkWCwFgsixQUFjEZYYQtcLViJasCWgCR5dqh18yEIDw+IV0Po
> IRW5mRyLIGxEx3UircgamdRRI1CRHndG2aDEoYsU4IBnFMSZZPU4Wn9pwH3ImdgiHDcbm3Fkazsy
> ZDgg7Sx07JPvEQ/iQBNhDR4OXOjfkrODW6McbsQ0oiMNcuoyc+N7OvifOgWku4C1HzRkl+CQehV4
> fJqT6rJ/BNtDt66yG1dPzysZTM7LMwbnEx+R2Ot4pz8v2Erh8oGungLDZsP1uwq85w5qNK+/GYxR
> lda+b5HNZ7BkTPwx8zQ0Z+ON75553LvCOjvNkDvo0axjWneLu7vC3d3d3cTnOnu4ur8iMYNQRtTL
> sbFCNVZyrYqMn7L+mN5oXH1EGbEAb0YaTxGHlqUiQoATYA+UUytttksRCDYVEF6Y3SbZBPNyBnp3
> W7SQqB7seId9hNqJ8ZE98zJFo9MRZaJDAWJLdBPEeI95TRX2d7VjUos4Tzlh7s3JDs8J3C7IFIGa
> YVhhgMMgFakF1MLSnEK0kG2pQ2yCLmy0NMWT6Cqq1MzCTEioshtbJa7J13JLBZtY1LUE8WvMnNb9
> AiQmIQUpmZhiVKmc6kgKRCUYmBNhNTs48oE04DB24K8G+WxWhkxrTEoYWhEwRUOlDSjOhIgInUbk
> lblTIYYU0DpBIREK2MNJqq+kHNMoXaWUvL00xW5gst+OwcFTFg1sKSbnO/Ix58m5066VT0VaQO1v
> 08a3UHJ4HI4PQO0bSZUHyJt2A9FGehe0aOqM1VUFDULqs9d48Hd3PUMMjE4Y1hrlIj9iQObhrYml
> wKbVjpIVhvubQkVFBtaoYwLSgXXL0B42GA6xBwmrArgsAnicAjRVyCJERAGDFD0chITJgYYko8kI
> shofaS+zQKmZWRqUFj0NkB4lU70fkH6Fd9AfOj9iAXf7SJsHdf7d+BNvO0aYxzOU6U65UsJjWvD2
> jqSkfryKLNrXR79j4viXPxGdxvodvUqgvAU41oSX3xHWohc74nHFuqjrNZbHNlsp4WpGGMkZTJkx
> 4NlMpSgKoqfMxwbbZxP2c541M78XbNJFasZmlIiZkmZtXvFpd6W7upu4rWIq7i6v8AmHNbOQTJJZ
> Qm1xtmm21tmOMQeBfuBU/EqKow4U+VNVQ1QHioXUOB+gA69ZJJJOOGepRtsQYQTTGIXp9R1Q378/
> XhmhnRsJ2ibsvdkqF1YQItDcGN/TGbVDWC0U1ZQ5WtdNJRa8fALUJCs9PdoW5XpB1VxHjfgGt4/p
> tZnrmmGYM0F+pWmMS2pIno5HO0yywLsgySSxyYLxSvyOTcXjR5XC7nBtyHG5SF+MYg1HInab7oEK
> T4dQZWHBb3kq7W1CL4IsZ8p0E4RJqhGhuFsyogOXBMXkHsVTntM8tsWxqp3CaHaFWLEwVczEvlK+
> BsZZcLbDJje0HxDPK20mo3hUWb4m0ZOOFVctYi6Arn29VuqcG54GC3PUx53FrKdEKLyIt52OHO0n
> hPBg9JyMh1o8kz35K1Qav8winemASMObN7we5RqIPtmRtueJEJhdPKhGc8lmMOZ1MIeKHJVtmPFy
> BQvImxheYBXgDXESKnCVpJOONvBxJwcGQzAV8BxVVFUWjhPiiEHhVntOpbg+0E2xHPK/BA2k0Vpt
> Zt0fz2Ug8RtscmuSoaDwqpkq8ihxqfUrtEQ6kA2EVpL13brrw1aWRV8xiGfGWkRxkeZaF8wQMNPs
> 24bgY9bPh1hx6Q05x3oUiCZDtZ0o6P0utlzmHJ3SMZyVGXI4Z1vkHfo2IVSWzn6t87a2nfe7SYah
> 3uYiDWol3u7stru7u7ubxjM1d1fvEYvJoiSzbYoStXidq2vbVbydH4GVdwztAHNl5lwjFRDAEzVF
> Y989aU0CUDIGMwyfIseobMS5ya7pKd7DFLoduAoWXSQCwiNBu7uJv5Z10sZy8rxI5GYSIxHbKzFP
> U9YscPmitIICM+wCHUGlkFmZ7GZfad+kPPgQ4XPrszMqC2shI4BwqIBEnrUaIBZ2EWWl++xGMNVI
> ssE1ghFEGEdNDIlgPQt1cMIlCYOV5IYV96sUSIUWc7vn0LgS7ySbd6s7853B3CeTSLY1p6Js5qyQ
> UyLxY10hpW+uRtOeGJUGMg2MbLCpRCfJ0nAHPr3OyeiRxNwHBmSS+8RKBhThbD1YgxgQLhsmwGJF
> FK+0PMSsEgXM8wvNmngTIDFerE2wtkMomhDKsKVVRMuF2qzKxloaznsqEuY6TYR9GKwp0bBadka8
> aNIpicCl6Ng+OhjcdOus4PIyeTcM0YF6KYoZ68HfkgUMnad9vFW+lycHtV8qqaKp3KHYb/d6OHhm
> XPXcE45nfz7/nh0Mz7N1jKtCqIDqqi5p6j8Jl94va6XKkKLYqzidlzlvjWmfMl5d9OSrqjawDvrh
> MqL8vdMYqsbZ333+PF2yMlDk2tstRKyTMystV3V3d3VXecvqKq7q/qEYwajbYkqhdixNYmYra81t
> wdxp3khXjktbv3M3nVV0wEyLvg4QmkMtq5zrhxk4pMKyJtquoSUNqoTrtoYpgNOCpppu4nbxLZ0E
> kTMQhbFIAwI2lifK7IpUGAqm77M4pRjmr5LD39J48CymzbFodiHjwO/gTKjrty3C8orIggnQ2w4y
> MFbaU2reBN14OxqWyEs0Ik2aDDsXgnHuHtlRWTmXl1aXkzIQZYwuE1pviCtydIl6EhMmFcdX44M8
> GSEg3LExqcVIzK81Nr7iWs/UWdd4r7ZmEamMRjCjyGGbJ6DR0thMVmD4jMZD26kSBg1Ww2kJMUJD
> 5GZcWkyAYynvnjWZ0zWFSNpvLIQ9+5G+VLghUUlIUC4J5SzEm8xLiUNTUyQzFxYFssoGxjwK2MpH
> JRuO84yS+2Paht5gbyGd9jz2d87mIC5KWOCCZOTYUOMuYOW2bWrGFKzkNCh+VEQN0ODsjflk6bZd
> aCDCchz7uAaCNiKwYxrsOzfJgHPMvsbjGR5fc8VgYhTSntQHuR+YmaYdrN1Z63yXeG0I8sarNudV
> j21xjSDjlXgiHgAxXObeDiD6Dt6zafDQtztxadukSJ0GThlab1zAQl3lkdpWuA9fW961vjHWd999
> /ddtlFwuKzK2O9mpFmdTE0YeipqqqZqritZqIu7v5iKYvbRsUbbbFiaxOtr2vGax2fXwp9LqrOx6
> n2feCUgDBIqnBdFRyVlcK3nqB7jr1Q5ct9j0GGojtU9UAmEOu/hWuFz0e0V1dJTQrBGd3K5RDjxX
> 1rAusaTf0k774F7E1o5HTRgPrHMDyHPezm+xKBM8+hQWvSW6MeW6fxK2ca65DroTXXglkrgc4Xc4
> NwxB4S08nh03gc7nkSa8+U8d2bnWVgIfGy8WUz2BiGpo0p5oM5jm2pjO/Ani2SamSJyWbwG2oXjl
> 9lCCSDxSWUmLcqF4q6whFsr1nAxlhjXj1vwVtoS0csfzUHZW24d6Lw0BB3d6LLZDi14hAXxHZ5Y5
> RPoXGejbxjc7V/rEWKBIYS5EeJHAoGxWIsJwhcFGKIoXG9ZwDv2XXRYbmjibw57OHIzr2eh67wXs
> /YrHWfCfYChExzLLszMyuDFaGhXaGKoZmRQzi4WpZXm+etNnc57cDm2DXILGwFOqP3uSAXbjodo1
> XV8OFuUiD/qh+n14azXaSJwue0jaC7oAsQTEMiZWp06fgzfhDeVLI/0HkfpLQJlD+ZPAEGAkrHSh
> WHdcb/uPrPXa+yipyZ2E4jrA555BYqnR+KsOFRURUQsp/zUhc8i9eyXT/9XuHFVVVDY/0cpJhpFk
> RixQXdwgfThrQoqqqmi8uDtGHEPTAA7YkJ0nsoAGECGeclFVJY0zBKKAG4Urr3eb0WF5L2E2+jCy
> L58qeWArnhSIF61OPRDAgnaLA9AthUaUjJoLxWCwAEXYLqLqFYi4FMckFLi2AUsLiLYUGC5iy2IC
> kFNVQc8xYLYED7D0aV2WLzfYw0nrdJ/Bh9phrR5vd816bKu6zki3i8ccLxuq7ReaTmJzGB6mbdCm
> umuj0a123+Fzql3ExJ0EDTDmivddnWqmrp04ysML3ng8PnPhH/kxbzOyz0lZAcjLrq+Mg8D4OZFV
> ISYuZV9XFOjIAchIV4EO9e2w5b7mfQ16cqHqSux7d3uYkOUxIiLHvhPeuCzY/OwL97Anhb3RP4Y1
> lsQOD4c4hzjpC2unX06T31dGDMVb3EPW2CUGL+XLkw9qZlxFDec/Rnw3qI8Xk9GB7nx593tez28r
> 5xhVi8uMJwie6ZxmPuwjNjEYPwYubNlcxNmZI/3kC/eI1msmt5GpGiQX/Qj3WNZSzv9kYxjH2evs
> qrwZ83Xly4ududpgDPR8eOPdbV393L1uQk5jmN4dYOEKY+Tq322WfLvhy1vYJl/eZuheXt2PbVpI
> oPP0Z7g5XEFCgCbCEhNre4XlmLz2cl8IfnJh1dY9dlcPm2EDzgB34PxY1NzWNUOCJ96/nddXDnOi
> e7W+rmwgi2AyNcH9jdULLorg0oP0dDjQuwZoWJpB559sBcj+jsqsVimzBAYSW/mi8vbLbh2ecSJE
> jHz8/OtBckhVdD373R9szNjPk5Njix3jFW5rwnx9/V+Hxe8bHt68DhDtxstn7HdCDaoc8O9vBRLY
> EC9jcYZGje/bjza5c+b82va2bqfI9IM2jKbCZ3TjKSRfPDsjaEh4eBTheuUONvXop8GCufsyv29r
> Wbz44RwYatVdsSQuTnn9vMrCkgc7e/v0sZ0MUTOMOmcYcx7RgleO5Rju64zjhk4jSU6W3ENTDQzS
> 8xhciGskN/HI/N+b56ZA/EYwI7WSrHs8yjETY7EdoUjnfRmFg3qmwcGBs14zK2Bg/BgfkArkQYxY
> REX2+78/nw/yw6g9WpAxAjEhJwqsCsiLEikio6QEN7EUM4gO2UJ0ekAMd7QQDHE5/PQB1dl5iKqq
> 9QQ4JDJ3D54H9cILggDhtiYFhaH0wQ2wMhTAoo3gFC1o5G+ag2gS+ODRiiLs0sq/ihBD/MoT+oMh
> EP6pVhDQzmm0EmWpBVESJGKRJJLtSQqB+9Fh/WEYuM/07tKCWIJDEQKUoEgMJQkVKISIkiLds0EQ
> irIKBeUQFfuRPeRIxjGMYxjGMYxjGMflrxkVLIkRIiRjGMYxjGMYxjGPZmODeEz1/OZk1pntKRk8
> eOIgrkGRSMiaGXHcoyqRkTQy47lGW4O1lbSUalHJhTHMpmpcTBolY4WUsqM1NFNOgVFVGGtU06SG
> JghRiMuQZUXFay0kSIQdyjLcgpYy0mS4xl6OIhoesoy0drLXEQg7lGWu1lTQyu5RlUIyV0rKXISt
> lUEUBMrx4xloO1lUIxJjGQo6V4xjKUdK8YxlLIllGyGKEGUUIyqEGUoOlYiDhBlG4x4xjIXIpSDK
> KEGQeNIYoQZCjpWTG4MZRQjKoQZSulZMY4MZShJRsQmMZR2UjJjHMYyqEGQrpWRQgyldKyNDhBlL
> klGyRwwpo1MxSxqYNErrKYxalNS4iDGkMg4iDGkMgOIhHEQg6VlK6UZS5hAiRBMQyGDiIJoY2kMh
> g4iEMHEQg4iCaGTAcRCGEmQhAlxjKpkGRjjZE0MuDuUZVIyJoYpkGRlyWMtcRBjuMplNAQ/1WH1M
> MKUjIIWwpGSWIsES2lsspSyWCUthApzgAlJFERFVghGAlhCpFgCJIosBYoKCqKlQRisQkRRkBGMR
> SKpIpGIMQRVUWCH/XLcwQhaACJsAhg/IP/YP1eebLoyDaLovA8H/Pi9gvOvfEcoOQ1nVWL4W0ta6
> ohICLl7asQEWq7jMV34pqQqOIimrUk+n35Mg7BcgQgYcWnGiOeWQGHua6nALFN3PQGy4lnDCGAf6
> 0WN4cKCUE8f5CwVqIqM26g8UePy+QENnIdiqoRcf3vvHmobfwhInw2a3iT1/WB+0P939z/sFJERT
> 8M4AXfnf0T3N7gEf+A/3gDk/n4A+irTjAKeTaFxbIh/EfxM7EB3bpoN1/i9grvkIAMUkAGQCH0xq
> GMBaAN8x5Ho3If7bQB3Y7DeCgicyf2CBRk9l7km9p6Usjk9F1bWFkYbZPrTW2Fe6yXTUxtxD/VhQ
> eVAYN40kFd4hE0Fgn/AghoTUv+yixL37y1rn6s2ziXHo2uFlhyoH7tqZ4a8lrJcaq96G6TEVjMbt
> kCtnjmEmtSYUxu14UEkk/n9InDA6Bof5IAfhkNgFOeZgShZRZm7AzED1nYDSQqGrikBUHGAj3Gbh
> BuYaUapZ8Ce0AuOvb8dBY5mHcfiHm/zqz4P3Q7giWL/jDvu5Jbng+Ae+yfpdHfGjck1zwHACrJq3
> 8+G8zm4XrDpInS2fLOnlKuF6eRiHPpB3KwzHRErpvY8YDhuMaKpfkWKIr7nsTL7ZEpRkQG51HOB7
> U910wHPOs2OmBhhfAgN2XvQ3Y4YFXwcCyPNgP19tcgIjmTqgdPBBUt47OFgfgOVH+vlNMANy8CVs
> 2PIsYseh9HL/jlQt15oeRO0ST2sTroHkleY0NKZSSlgk2AEVCGnh9hr/5c4favNQtEZIDH+z1+xV
> gqLCpsoeKW0GB9gHkLngnoPuhMTsQvEG+neeMJjvZ+Y3xxAmn5FyGBhIkA3kE9rPXCcNbPkq6QL1
> QvoMUCjqOlGAInKej4XyAHpPjT4HRJ7LJzCGQhPL7JWAsqqr4eE9tL6rYImCjMZmOEckWZd54/2c
> 79k1+5V/H5And8p8XnqZn6lIAvsgvNekoQ+ED0ELwG0Q1Ur5KBKyRPD4PmW/u5hlw4WGd/kO48Tu
> RRlh6ygF8Zrs2gajGAiAjAREiIwFWRViJBVRAVYCIkRGAqyKokEE2sPb+VOTOA+z7YTJ8udBxhRG
> bqk6kDR0ehGgu4u+YWyHkaA/ygbqRh5/5eu2/r1NmeJ6aB9byiKqqqiIooqqIwGagU5b0QWMk17L
> zFHG41MsDTOSRA5ScgzqBjKyFklTgJ+X2+p9y9LBggTqfJxTq7uBxKNZwWNznQU84XfKfRjSGYt3
> PoQvLek9od8jAsOClWvhfG3r8duhrWH2w9r0RKZN2bt8tjWpS2UtlLZS2UtlLSWMrLNWZbMjMZRK
> lLZS2UtlLZDyJ6YTvNtBqe2fZJHMCmAT3AG4FY6btm67zyKkdjCMiaaGlmmcyKkiiSaYGGYGEZoa
> Vl5zIqRIyMkMn6gdHeGpmqs0KkIUMIyGmhpZpnMipEisRN6HhDwV8Zrm7RRGA/49oYCBOLCwlkqL
> wqIlsfIN5Ts5LYQ5wNFRRhCIooo4s1qQ30Bv9smHP7fvR3OTD1cyqU6izts/4TACxiCTnIyqqqqq
> iiiiqqqidm/UBhOUhQ6uodG9pvNg/2LnKJ8hc+gR8QZwjY4IliHmF7XuS5whxB9ycJEHPkmJ9ou/
> VpO3sQO9to7vQ8D2zQ4lNqea9yjWMZ5jgnRyGuOra0txPb9GsRchj7GGb3p28Yl9BsBxI2u0BaWL
> NU9qBLPsqQogdEyTWPN3blotygHYmeBj3+F2rUmHqSYl7g7RwKgYVCDxBlS9spiMeGgZkMy5oljo
> LIBzq0VGQHtT2JgWQBXpuUm3OGQ2lk7idpNk5UX0DBjMC0rL25hBmkc1BeexVJBohdgGYdEYLcFl
> yVXfeCyPR9Q3veTds+Q8JkoGm8hfteOTFoMjIuEYkB/0CKAOyMCXCjyOebhLrJEBKxg+gupcBfYg
> B1DgnJOyn6lwRG9+ShsolQSiwav7js9H8lvw0YJ8/truYf1S/E4F2wzF/2YcbhuMuA0ocC4OBJlk
> ShrjQwNNqmxXaSTG0mazPARKk/f/g+CH8XuLr/3QDtpIvTp0h4c3iyIMjOnUq9Rbz0Jvz588OZuq
> pEVRUVSUkKGrIpqppzKFWlOkqLnS4CWBZw9xWW63i3ra3DYLeMk7YE4m1A4YLdlxWW6NXWnDQJEe
> 8DWt1rYVgt2WIa1ar3wtha9ijCmsUJtF3jyDRE4QXzEMKn9BsQkJC/JofiM9dH0t4RT7TnwdMNs9
> 35/7j8WimeCPojyRClKNCHL93X/w/n/1w/4z/ppZZY7vhhhOyuSkpKlKpGxsbbbW0tttRo20rSrb
> S3kcjkFqszbbNbZ5/I7SzebzeWWWcjkci1bau2Zjg8iMh4s7u7u7u+KVnr/hSboYgYoGwnf/yxwC
> AxgzYYRxrJSlkSKFBhk6v9aqqs9Ntt9EGqBmQQDgcDY5Bre/BKqSTU5WtVzA5tODxvbTQkFU9GS9
> m3vyzPAr7HBD8lHZwdnR4Oyjfgrjjwzu+trUp4s0/ThpG3sfVSXEhuDrcpENwunM1NfBevt89q20
> tL+4P4ifp1+X+r/bbX+kikYkikkJ0PgOjpGKQwwKSyRIlJaMFItMWMzbUH6opU2hsqAXRNyywRkY
> wcJoYscMSIQh98sCxb8JkCBZmZCyyhEAoBsG1n3gwNkD5CSSb5unGF1qLeMC7LSXRDWSBuDLwMCj
> AX2JxtwpPb4U6D5ALW+WpUlZP6jUaVC4JWFi4jAzKRAxq96D9AQAiJc5Y3bgJE9BNyEh7QJ/y8lR
> FVUViIq4TlzaKgE6h5jnOc6AVdcRZvK0L06HrbUMhcmGBKLAGDG5SgwxWbpgwswpqFlDho/r33Rb
> 5/yLVOJQOOQiSSBIiqqknQ2AJ6Z+Z+bFkUVjJYhBkBOIp/qJycx5AQx3yGbRR/QdE4Ox+80VtcS9
> LoLXQpYpULlIId3Hs5VbhIiggxYMGIQGMghDckaIjGRgsSEEgwihEgkid32zeABoHXQnM4Ce1pYg
> Uoh82E6gOV33LRpC1kmgDrBYQfxQrRQuqD0xCQkWd9k16sOrLnSFRLXsGafIkIT+vQDowB1PqOB4
> fAUm+7Hx6KGafb2bD9/Z90t/OfqIf43wvUvkQwL1YmsbKw5MOZIdEPMGBarGKqqqqoqI3+z9+HbH
> 3qgxBBvpd+NtfsCUosOuJgkggXwQcaojEVDILisFXJstQQYwIpA9tSCjkgfQeMkDIJgAaFq8cQHF
> kQBsMbtFIA7gkih9Yds3TQIKs5SAJDJAEChTbqEKULr2phmui5xpzcqIOdwoRT9hycSWgmBtaxdE
> Hj5ekKLDDuhU7/zaby9wHsQx/KhMDE/RoohQwKdyXSw3YHRqG4GC5k34boFokgsjCwUJRwcMnu6T
> /T3aHq8oarqA2mUDyfOLOztGOoLVEUpphTTGc0d1S0DQH0OiqoQxBiecTxK9H8PX8bp7Er9e3h+H
> n+cxzAdZrrazwSI1EqJUHzJqfU6/UFeZh3m2au4mutHZFOk2EjTvNuHjtW08cOZkegUHLoEKG5pe
> KaAEwsEZ56RZkRxv2EIiaH9WbG4XmJ/g2JMVyDaauuob7kN6ozxMR5EjA4lV40MSlwLyIcuMHftc
> Q5OTiAy0IxsEHG/gNgggU5N99R0cA/J3vsOc5BRg6IwDhviTwUHQavpa8EmcEn8aIiXO4prScMWx
> dkpwVxD35YTj9h/XjFm2+9LsYDHJuX4nUIr8Gw5gs4lRSjpl14z0Hp8iSXpXqel1YdY3HES7AlMW
> mY1B3KVmQuRcXAwclMnad4dDvPVynLi8sXMcYXEbbZWP92MaZ5QmX0dPN49Qh9nV1dV5g2FKkFLC
> 2JkH15kRxjDQ5zO82HihY8Y7YHORsbXX32tCEOYRIsmaI48iBK28s9APvBcy4E3sBBVpadLRzQtK
> iBWebIZCg+EZVAUJAQhEIxQHoPvV3YQBIX7fevkIF7y0qOo6ksM8cu7x6VdmVhaXane69EJN7EWk
> AoiQgLBf8vUy1NSjCJPYKLX9ck3p2bja8dySkOT1FX4TvK+FpwvJ2w2iUKFRYbmobt8Tbr5Ea5Rn
> pAPAqDc9IFCFoCuqYbgPPAUP/ggFQ+LyOQxQ4scBRxO5QvWyfl5HMZhTB10ucCbHsEdbAHG3X+0v
> AKHi2LBYMH/SfeMQSKQUikWZOBqUP5CP6QNT3jELGr2ZLnf+qQcV5RH4d1zWPGgAeMsQxEPMDNoI
> bEd9bdH7KLMOmkQDya94vvTH4Hfh2pco4JCFmwtlHwic4A7C5CIXA4KDmllAzUQfP5+fXli8AUt7
> jmnhOPdDQHzErV48miVBiRhCLVJIEoAjSFJEdcHUv2QADKBiQQX/eX5N7knw3PtwSs6BrlZDkT85
> KLinbsAi9Lhf8DwJAgC2b/NzAcMuZhxKmrbtV4fYG8ZSRmLFxMh2na4080BfuQCAl1VaYMCKmPOS
> EYSB1lDwy4TFK3qcVOcgUQjCBsMtD4TEEcAVQsdxl7EWHH3m3f4c9IFMCFK7Q/oQkiwIQgRrf0NR
> 4zQO+1WhO0KTPpoyedycIhllTmhjCpShqaHWC8RMz8J9CzkWHaxUd8QmVkJRsWywyMRC6AfzTI1K
> UaghAYCGLAKijhxvCkbCi96oc/nB4/QEfq4A3RGGgecHvlBDQ8ePivsZEkHelFHpaIkGhg0koexo
> osUWGOF0q7KCReRMUssJpwuPzfYfm+z8Cr7iXnXcfILUQwRVBr2E0SLuPBzEUoohwJRWDQUEMGNE
> x+wNEU8nvjdPtPoV8f00L9Du0NTuOh3H0y23KHA4Vw4WspI3gFQeY0E0NQunCkTEmmlLcYOhHjeI
> KuNz2nbRK657oWkAZj8+hcx0EI5wkYjn8JjpquBoXBnMoG8DiVLe3JhkLAwKEpY8UwxdZMyCiuxt
> LT7vpZyvCmWUCwwjiZZUsENUjZ9dipzPO5IujLZ9jbCbXmY+lDAzLsQsoTLwutC0oxmhwgXk/ejA
> 3NzG7EyyV28R8oTCFAmq65q+MnMXy7XwcOiG23JnaiShr8Nuc7lJlPOxyFF7p4SRshOTk6az4wPu
> Z3nre4kktnttqLGs5rmR2WZFu5llIG42kZnDjlx0XA4bVmfuASFxYXzCCOvwcfIeOEhCmEZJCYiK
> Uh7R4xiRKTbynU+pT1uOQrBvrMmOLqGXUvftiprYIKUCKbIY7Bw7LTicudfG5Paz3OPdF2gQiNyq
> JInIpCEGduQQTVBV1OXq+9zmdc7+uIZFhZMxT1MFC40DrUoNpdaUsVWOO2JZkWfXxAZgGTIZCbDm
> c6N8sZ64ji57A4cTU6B6AiHcsSLA2jmL2pIkN9IiNi5jbRYYuccg5YZoSLGEbqhz7DzIpxKPf+g4
> P1v0Ip0Khs5aFQlBSmFNAfSq3QHe7Z5Y/4EDvfrgUEDQHxsEBYWL0FPvwO0X1C48DUMAfEXcNSw8
> TSD4ROYL9ns8a1+0EXkgIyAg+Z6ieKFaaGa2Smg3qB+EvqHD2GXlO3AwhNE+VfhyXffXgOKTWSqn
> Apsbc3LeGS+HASgNRw7QQeta009T3LYJEJHvKzBQmdnAHDmRe2Ndp9UCZtGBzn7jOv5WvGQ6OjfC
> a657FJDbRRQx3RsqYGNxcWRel95ILxfgsDAvMTiVluuUTwzNqy8yoeOecbeN1PPDFI3ext66JUO9
> hzBHcX3EG4T6Fwh5PVsjTeMnJAOcbIqKUcEYpsEHoa5DY3LksChnWXUGxxGl1CsFtmFkPr64nITR
> rYWQy0LnDA7d5zM+Ce7b3H8gQ/HYbJPcsMcSANWboue+Bv16EOvZtTYq1UTrUnVRRJBxKmRiye86
> 01xUVqRNMyKhjgNd/RToObl1kmWJYFCrn0DORg/jj8RjIqQKaoBgHPXrXQgZWOnC7kVIFWISHEgd
> CHAAV4UxC+819zjnZhqYu05kPCxh+1njJ6aKO/8hMQ/1blJEIGWLYsIfHR7xYD5J3Iwue1aD3pO+
> 1o3IkFn/XhkkMISHt9/IRAX4JWMSSMRFSEFCBCQjFPHSqk3dvYbiukxudR1eTtxMA+oOfLt6i8j6
> wYnpm0ZFkCV8rB2iw3AB2R87E0IZ23KFK5umg6MCOA4zfF3s7MikAUcBaH6/noJGJ4OUU+IRE5+G
> nBVT9JAM/mq7UvAvcCx1IEh96wiv0oH3kmodu2g+kQMD+0KAg9kIiO7DyssnWyc/Q20pSFq0djYS
> wK7QTuAE6/GyHrYgw+UAp2AeeFNCWbpmBZgMswRhVSlCydkkkOyLBARg0ZolKFB7D91BwxdJwEtD
> GFiHhAns9BhYk0BXWVRSxitEaIBXcbgE6RAIkUc6+F6WYvE9DsPk0r/BV8adExIwSC1ESlCgYoQQ
> J5QKEdowIMCAhIJECKEIAsbqF7KVAFJEocNsC7iNuQHN/i/y28BLtDmSCz5PRlBmYIpEiG/EK8/v
> OoBt2zvnx2q1eZsqVznrqmSFiqcSYZEcm7gqhoskPsQw2GMG6skv+52SZqnrRcfC14TqSGAKcDKi
> 6sShWLmBcMLYTi1ByhKMOW0gXCMhu+Xy2+G757xNxu9JUiJJIgdEgG1d4YeLBJgwZDsTZwcXzUD5
> 7hjgs+NCDkMnCvQMTElYfCxH0D1EnI232VbwvScvfx92h7bDMzKBIayd4aqWs0V9DFgaBME1paPk
> Z5hqIsLW8cTFDoK7JDoXZ6MFY7Ojo7eU7uDc2Dg+QiKPe42BMZRnffKSZWGJDA2WTEcCSopJ7LRr
> syMMyMVo4QY4oEh8J22GZE1stiWlCtFxxEgXZ0ObOdm4oTQhkKRSzboRK2wHYfmRPgl6JxG3XTiu
> aKhyZ8fEZamaPXBXsEFOKQE5xRWwPQiBwFYKlxaAAoRpiIdIBkfDn4v+MAWazJmvQR1h1LQuJaSJ
> VGxaudxUFiM80IUxdntYf5wZn/j9dFA4cHbh2kBIv44HcLM3W52B2E/wEV1HcWqwpWc+vn1106hK
> nX18b9AU6NyrPIp5/QooYdMHPgwkJJKAAdkUOCgbSxlJeiAF4ic73ZukVwL9H0O/gDk1XNZdxucu
> ZTUtI1BGALiGxJTEgGTIRsv6QCsZgt19+xE9Ux0CbDgvBtQ95eC54RuJotVjlbZiqDlzOddTmVYv
> 14BR+lVMP4RAYQHiBoVoRQioHGPtAUpG570iH4heiOPomew12XNfRqOqhHlFZiYmf7QA8ApLAKHn
> UNcBwQHRA+zv3jrELCym0UPMWQ1onZAMYAP2ir/GlLa222lpaq21bL7aFyltgVa2aAJeROsnWSEn
> mAPuM6pAJAL8ZAInwcNxKRmmVPikmuYF0NQXqgNQX5+R19Xw/DIfA3p1uezfzz+vnv0DOS8KmMdO
> VDECMxr2GT84HqJb/UrTCCbCvgtCqDFBDxGAFzy9mv76q8mRW5N2xRvI2BjIYpcdACr9mBzvS4rA
> gzIOBSMyz9G0aMTb0npuV+uz88cI4ooyLn8jHHOtrb9eilFXm6ZuvjXlxFPEjSPEVrPC6jkMVRm0
> HVQVUAEsYt0xELRplPGkRo7aefDSU0gIgNpzD4w+2yCg5pdqMJOlNmXVIqWqq8MBUFit562ern0V
> OLoMh18uRo3HPT14YMFI9VOQWbXKSFJMiM0DmCfSh3Bst2BrYKpNMaXKVMjrJlVFF0GF4vjG+QiG
> NKb9ZwQo2RRN12t9LSImtDKooLCif7esZpaajb6FrQ+R8/AULpz48h3lzyavRTYi3lRI3bDTgXCh
> LVVtaQuvEihtM4I1yWkWJVmL5jj2DRQeKG+vfSq+Hwdpb8J9Qw/RJdG2HQKwrBJULugzDZLrwJw1
> sdrnpJe6AKIKEmQ1FFBZPNYDYssYAWYgJwvr2O1QeE8Qa8+Hsw4P+D07fecjlQ/IPunXP21TsaKl
> aqqxjUokqiRhKhCSM7DuPlt6uf5vf1+Pj6Pv569/y3T8WvFxz4eDbiaHlMw1KPp2vke5aJGSSSSM
> IREFFFVRRfHpyOvoeXQ6zrPb5ezYPKcTf64dRmxDund3RncFSCCRBUVRRVWZmZmZiJZl0Ncolgab
> rQwkDBp2D2EDY4EYqnZ7hfch4sVfbds9plRmauvXdI+3fJsjba9IRg7s7KqXSqcl5MKFUn1jSt9j
> MxOL7aA+6SqOO0ypEyrfiGImXGxUSLAoOgE2SSVtXF7AXc8FDdiHJXJya5RDHUn57yYDxPR2NsGe
> Ak7fwOdHbh1SaPA8Jlk56ncx3WNpOyshZ0az8nO6nacsC0McG4WIczJ3WsJYhVGLtaHzwSIxBixk
> CMQYxiMYMYxgJGMYkFABQIjHQBJ81Q3gsgp6AA9CAHJiRRVUVCVQjOYVhyyduQO/TYTvisbwbDcy
> 04RchFKeJ59xUV1tIsgixB1KseYwYuFBFzRhq1VQpCw4Kg2TzmiXUSFlIFCrQkBJBCgDHNAbzxUD
> Aigkk4hIzY5dYOcmRxCPLiT8vVi5GP7BxT2HSW3wWgceILnCc2AJ6aqxc9odvcYTPt8Ch3h4eJyN
> gVy/SdYPZh06GF1+ihjfl3FV6eAMr/NIOYuvoAfnHAt1r9gxeogyIpQEZ4qBbJ5CpBLkUIqEVtE7
> kIkAf+IjQUNkgcipO3ubxZY4OqF0IEwXYzQpV4zXsGyrHq7rbPl9u0xPGA2cwKinag28M79MjwCK
> hnz9ucp4OoWKwfBV1S6juPQ6YB5nT0mcDuITQsRJWguWbgHMiJJAS3FCwwS6JyC+94dpkoppbtw9
> 7aoMIBy7yhps5xVAfiLZmHPGbg2G4xKTu8qXE40B7jcHRD6/PEdiXOnhf06Q14DrJhkc0B2ys7Or
> b4/Hvf3Xdj8bxLwEunppx+Yi/rQAsdCw3F1xIkJLoOdhuIy6xJAvIj2F08RJIDy+trgR+5iJr0fM
> FH0G6piEfuHEgnmkvMH8hLqEg9feIQeQu0UkFzIMe/Lz8u3yePn469X6/AS8a++3yjzlF7osUEu/
> qY4eMZdkLK5n8wlQEuIgOxeH/Hpe0l/gO/9/lqRogQgCSJACKeUWARiSKMZz+P3v6j9RtVQrUPs9
> P2cT8tDwxQvgjk+snLbGTZnTexTJcseotgFLDXRFQOEdPdP8vbbh/NlBEPuPxCnV23eghGIyIxET
> CJihBn9E/PBjGQP80VBYjBEF+JDiUPT4UeXcLE+UOiJxEvArlFogItmrNySiQxOqKkIMFYIIbOEW
> cCIsoh+E9bCH4HbFwBvw84iWW2nswFBw+lU5ucMg5n2iepoT/5hooHj3gOdOBv9KKevUw99P1T5p
> YHwBNoh1O6GgtaFEKd45Dq+KWRCnLEEDpZC7a3wmwRsfhv0CxowmtJWmVjzT4RNf1WFxgWggmuI0
> cvNRxgGUBYRJEQIEBI/hyFsQ9aKgSi5/Du6/syxHYcaqjS5sOD++LIBBIQQkfbYUy15BndVDdAQS
> g3qwDWAAFRCCJwgUgG+CXq4EBSGlGuykNvziOXLePK/MxtZ7gkRjrfaRpC+qxkxGsayxS/Q134rJ
> oyppCRcq+3ljlY4YdtdYXNM7JhUbfKOSqKNisC7GSIsFbEC7cUxt3TS+uMvbOmkIwThKCyrFGCSI
> Q11JgBaaLJvxIZMZGIhgb22lKJRWFq0hHFQwA13/UWF0PEDgBrzoyQhGJCSIkgpFiSQGQARFVYow
> RWLBiEDNzCkbhw0osLpjyJ/CEGmqIqEgg3FyxUPs8AcQ6wryAWbQjImJCohUmvZ3gaDY2YEQElpI
> RCkTXB4OzDoUHZsN8g4rygCyBIAdrMEdXWKavJU6RkY19wnqVT6w5m3KtT69xqKode4TqOJActeP
> eqNlFYRIRhCIwIIAMCgqhWiBoDrnBJGE91iijyGQXy9aTFD7cxYmgAjiaQ5zI9w+DgF8j7iPkfqH
> Br2pmUiYTiPYb8JfJVkCRQLAZ/fgbsVYCyAlWwKIZIDGoCSJRAK1+3r5fcb+GPICRDlTRr70/tO0
> +Y6vNWv3QaR61VVVV/GwFoH+JQta3qOP9ndA+n3vwh+A47z6/6gLk+clGPzno0Tjfo36Jz+rOTfg
> OVyfPxe2xprMwZFsx4TBucngYYdOu27N+s4N62RxaXkNDnRvg5Eo7I0cHRzRwc8DPIa37wO0yXvA
> UnHjw/gsLcnoJDk43NdcpxxyU788JyaQbO3Bzc7mEU58YM8JyYkN8zo8Ewxw0GDbRgUAnURBZUT2
> lSzJWGRjgGKFzJIVgkgURD6vdtsbHIzOZx58TBuHE7XyxkQ+BFkiSCHv9Es0IiIjGB5/jQoP6DJ+
> OWeH6UPfev2zzw1lxchTQQs0pxHQYeIMHr50U6zrz/Qp3d4ejAPDzjkW8iayBIPn8FodvFd9UWXT
> sDYqH7eiSDAyDI4d87bTYwCki0RiRMLbTCotNt7fOq79tMLBmApzcp6FF0sJtr1aiqB7stf0GcKG
> ZptiNlDINg/aINBAdOPP4AJFDVeZ/DjR4FUxOFghygnKfggDHzCIRBX08lQKU1UMWBQB60lhTggH
> lDqFwvMl9G+N0RmuAd8Xz+b8bWkYwjCJhVAmEQKYkMJXuU4ufyALmcF8ouhcQtrHQD37CIfMDEcB
> 4E6wkE9ZameRy4+WbC4ifLJtPZR+AibFMgfAGg8hz+glK3AiStApdkgUYExAhoPM4GBPmB+YBSbe
> qsTntR8KA6I5IDgIhSqYor9bBQU1OHigdcuntFCg+ioaLvNlfWGsHHQ4xJELIJFQb3grcwA4cngB
> 6JU1vlqQgOg4EnAUFQSx252QbDESg1fmdvc2ChGoDvgJoMSAhAreEVQbyjslnAuDgBzVnxFMgA+o
> OP3h11CUCPiaIlwTgUgeKoYgcBC0eaUrinS7FPiVD8FQ3eLqPLQMgQ8vCtNwYXoUWwoJaf60MVRR
> JvaWgX2EDbBknzIyaH2/+Kr/EURBGiucSG5ZuSIGREDziQDYfOpQ6B6hTqH5VQ1Vw3SQgyOoiFRg
> BpEFEntpFB8DsMBkUmNkEK9OSrH7XKwXoUql+menkNP3wP1mClxYhulE7cDeiDIqkihrAKghIkYr
> IjIXXFbndspgwhAmwrMfqIMZECRQikFgjASRBFGIoyIoxFGQRCRYCMgAskkCKAyJCAqxCJ1hUAgR
> VB3fuFxXGycHb9dg4bdUxg3PZhuAQsYQoNwgWRMTS98IZ1BhJEgxhkYc4/sEQxwg9SqVrEeMgvBv
> Jr6lDoMDZlRrOMeoFc/R6fa8Y/oHk7PrMzcs95UVg9gQLEieXEx7vsw/JVBKoTIDM/WVf3Q6gwMb
> 9hdX1oeCuq0mBkwPnDVCET4/k0AoDhR+U6S94sTsYog0F69d/yIjYcNC1BI6kcLlxoMMkSyfxBxG
> 0kCQirCKBGEEOSJHCfMnnF7PAe0hv5mHtRT1i9QorwmO0OP4D3xrIXzeY5CEOwzYHb0JDjBOsQCP
> Dqfd9KrBH+UWWtp6kPiwyKojwBIchhOuH7h0NeBfiAeyBPVdGMhP2yoQQkRiHeBlliilJ8njQ0BB
> 65uiiEgIsLhRSJA4+BJMfgpzFeQXRc9gQ0t04v3AQUphVHKFA2VVDtDPCcBJRR6RppmUahZkGJeP
> hcLbrpeWEujyWYoFFJQWQYcy0A0SVskGSCwRPFKiFUUJmNJc/jyRy3fg+VVTzMbGYmU0+CoNwxhk
> L8K+EEAH3z7GJEWEkEDcnoYpq8Usqbm74gIxQ5Mx9HNgaSahhuw7ZFVfittW2rZaq2yqqqtpXR9f
> s51tzATt4sJOd/ft1Tl4sMzVVDUjCcSqk6TA3oY2zbZ2JbLyYaHo8d8cjNXe0mF7gqhYUGBAeHYP
> YVzb5RzkmxDWKFiSSeUnXHdeofYswI4E3YXroQuXzIZxbY0qrDmApVqSMR0qghKoVEp0wRxJ5OoG
> pMXW9uqsA/v0YyfUE2CFAxKDOZQnmGHAKMEgmGGwJvoMqAe9Nhf0C8Be+73kIQiRBGSetliCkFUl
> GSpAiSMAgsFYiLEBILARqCm9NQzSBR5sAf62PRufZiLvIO9u7QDpeoGJIj1/vBfGSdggfCCqMFVk
> RH3iURVFCHzoda93D4+V0OAUU4xAoNfMmMKjTtoWhjCSDCAkv/c1sMDTIP5QqO3Eqy+lMVDxpcje
> ggNnPm6UHAQ8j2skh3RklwH64o9qGJAZOZp6lGdQ8wIfbLegf3UWygtJzM2IdK8sRFsV22HK929w
> Ru+6hQhX5lm2Al/JSGuRmgjrnbi54c7jGA+UFYEFgAe1FPAptUPm+tV5R98QY7BEMRPsET8ADUQN
> hxIeb8UbPTzv6+hyJ5pAkQkIYDIIwaAUdGk8oIXFNdFD4/lpzFi64fUFkLRGoR2urzY45G4vHeoU
> cW/OPaoX4T5AD7UB2LuV9inBmqD4aaunsU4DvETDAnc/aKRAYgFetc1DuXzlv1VnTCYfuQxZEDSi
> nZQTHAxGMQYiCHyso8EzoUA1qEnZ/OnqMzMTBXI5UcVbQ1rqjZLBIETVGgF3KoWFoF/sRgLhDhBZ
> LxRoIpGFURYKbAUBdVYr5YqUjEjEhGAQdFEE3wOnfA+1IixDmVUfbzKGDgPDGBsgWdKjYgtuXrRS
> grz7OXh3AHGu4J2TYsE6A6owjCPGmoxQvLRwhFR4AZRPM94DMB3AhqpGrc5QzemK5JMIZYbm8nGq
> CbBFJQpJYG0LZgyBvcMoTcZJohA8e6CKkWRYqkm+geCpUlOfKBoFOLKJyhWc5DQIthLCksFbLjrx
> TLjMQ8mAYnXUCQGkTUmuZTKBMs1gG6UfydiKnGqBF9AiEfOgPaG8DMOHCjHYM493XQDK/0ju3Qq1
> /Z/LmtHV3s7hZaNlLZfVmZC0lAQkPAWIJQEAIhNWihXAWUNQhA5ssEg7uxcMwv+85ZfZqsUXIxqs
> SZ9EFYwHSLKLj+QHGKUpcYQCQiR0KqQTrRxEsvR7bB/MYi4YAFUaoIxVXqJuQRLbXHP4vFT7Mtxh
> 21ejyWAx+WSWCBUfiw5KRUAHrP3CJ4bGd8SFnsKalMgRCoQCyYAF7MGA1auzh6bAJ41XaImSPpQH
> QrQxSAllbAlgo1Nnw+Ou3ooL4UEYpFgQkBPzr+d9RTAgWEOQo6iKcA0Q7XSDCAdrsoUsPJA8H1DE
> li2w2BQAduGEgwnMglkQpEpAYBQiWhSyRJKAliUiUgwKQSgykSwQoCWCFBKNgliWCUaCWtwuoBjR
> QgUBy5gPJR1VSioj2KxKUC21sSsgx9sKWAIEkQNx6QbKZ2B0YMmDsHRGj7OkgVhm2LGxqzs2xS0N
> v56fdo2gbBLAUaWv6fXZykQEE+5mIGPFCrRL/R/Zge0iHQqEOlmHoCcF7G5HpF2i+C+kkpH9iA3F
> tcnQFrFpY/8WuiHT0i17U9Ym/vIXih0ob2lgJmnHwzsuad7pQszAeMBSChSG9dhCAB/CF0JdlQbF
> hsoh4vsEx7fBi31UMd7VVBEN58wl6BydmXyuCBlKIIblAtFZaB1TITBIib2gG2qYI7FOtmoMBE0J
> oijS2Jma42hix0lN8ijvoC6QTC+BYuxkX04WwuhCJTjKUIkbaYKZcUEF3OlJjNOwyOFP/0EmJSEG
> GAIyKREgCRFCJJARcEs4KWQBVKMCvxSQpIkiSIiAMiMWSIvJJCokSRhGfjADH5B+BAc1fyDEcoqE
> CEkkiQcWmcabDGbikBqAiSKJb1nrPwttiGJgDskGJgtxwAwFFaGCCQIoSChIThKArVcjLJEWz3lS
> IYI5HEnHbFFQN3AWhsojPEwAC4WRIJUCgjNEEmzMEVQNaote+IaYMgTWwKoNImAsYiNxbV6xgNAA
> jywUDBxGCoYiIZAjAC6MSAXlMGpLJCwbICwHoESYIHfbIHWFFuCMyAv741SPr69oZkC4ULqlZqmo
> pssqLYkkIIsIiCN1A0t8bYBQuyRVP50DgIWa2i2sZXCuCjguG0GxhZD/WRSxWzyXsGHehiU+8gPD
> +cB5hy9YPURkCQP1n6s/FZA/C6FiFPppS0taHrXsaKaaA9ELEGiEQgBxGxE5mqRkEBfsiA96Nugg
> VZgtCRKAtBQk7gQJHAwCXuKSYqIMkUJNrWJRGoqGBFEKCKWiUrgNKJUUaYwkggVAEYCLEAjcgkY0
> kwiELNoJIUQokAKiMAkyrIKkIMIMZCKwUijEgjgGCS6+CkAjAVYiQALpBEuNUBqkokBkgRh2nRgj
> D0Fge9BEKmADgImQ0BmB1ghw48ROPbCQuWIkfck2HWyKbHmBmobpgAm20UobR3GTYOzt7/uP8h7/
> EFKg3aMMyD97B/K7h5sRn2naQEaIbbWKCwwcSdj5tjnNjsnw9Vh51YK9oI5e4+X0HVbogZkAXcQE
> eIkUFvj9XnpddYMiz2qH7VDQTpvF5w0AG0eCbXXueG7uKK503nEeGSHp/cjW8gZty3GQJJJJ/Psa
> qA8Aj/ZQDY7leAFPIOQcASyfainl+/6wIEFiyEn20kKpEIqRggxGMRZCKoEChWMgDiIPnTEDx5Kg
> UB0sYpgXf7xBCMkUIkGAMkCZzggiiArBAggyIkIgkYCMYCjGBESSRGQBBkFgoIisYhIgSIEEUEFg
> kgDEYgqpGDGRBJBGAAiAgZvRFdxqoOqj417+Xkt0o4Av8wIv4gxFIEWEZBkGEQSSQkDU7VQEGRhE
> QSByKUIWyeDkLTDApBIwWUYVgCpRawYHY5zBqk8Wg+2wgWKJJBRSKipILEh1F1MkLEgIyIqKoJId
> ZSRjUEBFSQYICwEIQEJEYsIwFSAEQEkFDpVpRaAihCIsXe1Ak0HmEOXAcH+AuGABxCPlTcTKoQ5o
> xNgacrA3yqpAPOwAqpCGmpbVFgiKSIQYQSDICMUBlKFxBXAWAJLlGrELEgTPxDqfuUAIoWRISKjm
> B1GGg34jJyQIrtRzKlPq7CR6bwR1aAogB3FKJESDIoznQsU7QZswDTCn9Ps31FkScXZwmUhIBSBE
> dyAcYNKpmQRrMCgBZJH9RoCrMEKJsTZ80Hsuoe/8zoORp8D/KMrBB2IB6AiQiQgKRiNAjtB53nDj
> HrXiDbWgN1EN1nhTIgUh/UIFgOU5S1Cmjxo8IOAwvR8Rq333ecf6YB2qK+VOqxyUHgLK9KAfxAiO
> 1C/6HYkCCgeg+nwHjA4hAnI7ih1T1hn4ZQ/RnyGtVF+cP8tA+IXulAqm9D1lNow2roWj6jzVHTQE
> I0caGfNPnuZH9ppSBbI0KS7djG7wdvlOA8J5PWVHqL836PkCB5zu8LMIX0UAp4EPWw+h/mK0kTOj
> IclijgccjIARuj60Oa6dDykeUXMRDDBPqGqqEMSiiV4jHRvHZVBmQQuRhUKHhSIn5p2q4nYc6hkL
> ZsLARYKDQVwRORh9v3YAKVrDfCE1GAAUopzoDbYWW64eGtKKqxgput1wsCZeWGSofSgOnvG8yJvG
> gf/TF4v5rF9uvGgwBXADvAsjSd7A+4aULKpaC4IwWn6kERS1/UkT2kB7vfQRyiBchBAdhEPC0BIB
> IBaPjCqIB2gU1lHRM027McEqVYqqqqqqqIqqyCvUwpwSWSk/pGVjIs5sCF+KeHZRgGrjZCQcgCUg
> wkigJFdhSgFRFmuSWfrzWCRNi0RVYqqqqIqqAVEpkghOgcwOtwC/M/ECPliqewRDynzeQMBAP2D5
> e64H61WYAK4EcNCnfUNUh0F7CD6Ko6YdYAUSQe3/7RURzVCFTzC81G9ykCf2PwA2ttNETIfwVcw2
> M6DrAUfUekMlEOoCb2AEV5AsJzP4qp/cr2nyopcNB+NRMH+IidgPlQ63aJoKQ5QDrPeCwQgHekXx
> /K18Ya2tEU31JyioIowFGCQD3lkoyDCMYIogsnvZVQWIDEZEikWIArBFjCONGfGyS6sowHTJUUGC
> qxCKrEUVIPkZgojJDIeyMmMLBBgS1OjZTqYryI/oPhVTFAPWBg0DGMiRgxgZhAgBEdCKWFDzMB+6
> wt9FC5C3KFyFLlNakRIkgRYMxXZVWVBgqHuY/W/fG6lKh1pKPoFwBSxMjwEUtFIrBCQBHAHIUXzH
> MRO7XNASA5clazVjeqzAJQwc2LlTIsIluCTDDmSBhCFjCJBB2Tu+ohW02RnxX0wYcnZB3lVqKc6h
> 8Tyapw7z/mi9iA0bhPUgF0W4BcEuFkWyA2LCTgacN1w/pJgioGAgBIKDzH9EAvuQ3bY22pDz74fl
> pdA7z+gUc4P2f03sJ54feL8haEANgOgc/SIhOqWgCHbAK+5OdVEvJAkjADsg34/EZhbP4C6hLodo
> d5AgjtlvsqUEzrQv6bjfM/KZPPlityAGJSRVRpoQ+lGKAlhNUAhQMRGnUyEzc4bxIRgkUhFQ4PwI
> LMEQHrCIgpsKBpEFuZkwoMJdAs9ebbM7amqS/g9hWt0JGwkfhPgoECB3UFrKu1i9hVuEAFOggXKs
> 0n8loLQooaYAKWCAUsQC6wEoCWiozyAocl6IEgAQvcH4CkMJuoackDASRKKONIKGBGQSSoLTUQFI
> kVCcwIe0wtMfVzViGJwQx7nbJrLiUg3mmc+mZNgTnyOyGHCHURYmT7wpwEn9FYbGjRoF6Fx3G1cM
> WsWfnlN6OItdy6EM3+dx4YUQlnZjEi94P0LuDbDyXqFqZvyC7m7G7I2rdbVJVB15V3IMtgxy8i40
> hTlotbUZVKNZg+vpTcduXazjsdk4/c9rgXKptVj8KjCq7GlNWps6vgddGzky6HejkzIweq2RptAw
> SGTJZkJMBSoYd0rAwsNzmZvyCTbguENJCaq2rILEWXxDUsmEuYY2egIIWAgUKGSeAQTaOAXHEzxz
> VNaDwtxJlTvWbo+yeYtj4nK4kqaFWRMOgdAzZRAZ3wpwwMDgADnAw08AYEzBnr2Z2smkBELlGJjg
> xqXU69rElrrm7tBmsRx0YyGSQhkkJkkG2NsbY2xTDBDEQRBsUUBBDj0oCyRFuoQDFyVC+CuKbluL
> +G+3H8wuyrWnpTVHVvp3eL6CEhCKoqSIEAbEI28546fGfzPLuGMv5DzQJf6gwAkFFYEwWpbbUQFz
> sR7qRr1KHc7ufUUepRVQgd4/gpcU0ggHZABDQBRg+XeyAKCkRgKAoCiyAsIqkXoqDrzFgTIu3jiC
> hazyFnh2+YOjG+FU0VmR+YEIBH71VQjgqGwuhSCCcC55AS5HpQLUEAN9mhSC36hotEWRQTogOAge
> LuNnxgEgSEgkgyAoCgpFFixR5kLCB8SHQ5hPbDQDeHTNPoXe0B60e77rWLWtYtalDzsjIzzvhCJs
> xVP0/kv/amwcVQ94+2ecJNmKRQlKhBaPeGdsFFTbtDCBKTqbAqVKEdF4mgwkmk2KhsXQhhGDEipk
> 2BChhSyUZIIu4ghB2WaKBjAwSCESWSm6uDRSAsFGD4dOYO8PyJhP3BlGCd5jcWDGYKrmUFUQcLhm
> Mo4Zhcctgqk9TFiKqoLFWIxYMGKCwBH+QZL5+rzn1duJAwOxUG1jzPNHdtri84wUYacdxYeChHyE
> 38abQ7gTuEQ/ah92OfHv6B7i557XD4aDbPHC9zOylkyRw8C5Zhl1LY4l8hUxCKiOwKryMaU6KaRY
> MsA4CBaggqMUti4REtUZZb66p5UG1KmTLBS6rBnCCJsqjAMf31NJy2EdssgM9GI0gaUKi1bdsyqA
> 1Cpr8CYvfLNpsflUQcUvSb7YdCxwc1tlgfdyR5IYhGdmxIqHhNLVZz2JOs4KkJeDkwwI6YrTvaal
> 2woJ2n8dYM8bg7ogbCCiAymcVMdCFCwgiRQmENDNYEoyINwruspHPe4MqV5gyHC6dL0Rrk9w2Y8w
> NbJ3l5EmXa4SWSjEZEElOG+UQshLbWJP2q/mfgD2ZbYu8gJrarmtcbVL2LNECspNZMO2aORI4EWA
> s3KCGi+hIGpDQpMBoFoYRg2QQIOUtiQWkHiUIi4FmrrQBcisIqySDcEWCpdiREvSgS7cVgpHqgDk
> xAG477u4d+6AHUILc1OZqCeJZ/GJPGdv46xzLjmWQPKLFjwTn0pDz4IsQRHqnIDZxDFUUgODIpIQ
> tGok50GMsVGmWpt397ZtHMWbnJMyKMa6dblsCrDoQqxWNJ3Xq+MVPfpZUJfY9TJ0pOWVG0OkP2Jv
> cScolBDzszNLWKqnZYhaBeePZRkRduAtBYzNHG4reCsRhehkECRURSMiwYCIwZFBBIUOwNd3InDy
> mKCAa8IBCQIqaOBqSff3Vaw2GMD6E/XMm2betUVs0ZCyAeGG4mENAMQJBKITKWHcbgBShAHUXUE2
> A2+xIqdsozIpkKDsgSEikgrIOyJz8zawXvS1H4Zl04B7g6o8cg12BguuQGkXhDlwEsglyEGMUIAa
> XM/MC+K+BGdOzgp6m9IvVBhhuJZFGCIIyJ1WrBrS0oyKBy1RyQsggwYgk0liqUSwhEUisVihMCqR
> MRUGgqwIY+0RDC81DnGo1BDsg1LRhAMzBK4AB3PrDMSPAYadr6DUA1IC8kkiwsQFIIJBBikjIyMA
> Ug01QTIBGhhFDUWwMBegdDoGL0RLoKaAHCIbBOTFHoTHuAUOK2SJ/0ulKZl+k2HdfjaygFvOkdAL
> wqNLrB7wfr1e+nIzf17bOt1IgqyxDQcBpIiYjduMgKYqSz0uKDRdihbEvkZminZFC0UBjmoHCu42
> G88pFJEVYRigSXwfO/JGd5VPlmXDApSzBDJhQtopEtEsoVsGAiUs9JhdH8RIf5nyw+YGV6UB2Kve
> j6Q9KRISJBSSCkRHxsThOX7gjIsGQInaAAZgeE0+dHoAtFRkD8e5BGh4vi1CFkRbayhNh3BtfpIc
> AbtzkJ/jzI+5cnpPQTzdZYghp5jjOEEPj51EE8O5gu0OZfzggw8+SB0qHOaNMpDwov8IsjzguEFC
> WpUG005ZofrQIOqoFhHxI/f7MsEM8H061L8yICcCJhEG8UQqFPypayrcjgUp0Za50GoCj6DCAfVD
> F96odQ9oaGZoKDnLMLBA7iynEN0yYKv9Fi9wOvZbR1GJp7VMQeI+g9/pD/K4P1kjDyE7g0ElHBQU
> RgxiZAn10AKH2QBQ+QT3I7QP7keYG6qe0cTgYcDDihISEwGgACT4gnVD5X6UA5ucgWkCT80iAaGW
> R84epGghuShz3zN5hE/x3wuwA7ACCwGADFGKKZId5ilBoiLF3kglUVsLVsgeoSfQgQiDDTy8Ca/3
> Kh9x2OUnojmGWQBxgAl1DQ06J3CfqQ8i7ID9o4ZhR+cI0knmD0AhT1+N3+dDBIJH9PBLhGAwWCyI
> DhCiCfufsQG2O+j9KKe0RMR5dhk6ngoLiSHOahLFETYRHzmxAOcEPkV5lDB1+Jo9kgO2uiG6fQ+P
> i9l+dIRAj4lyoMDkaO2K1HrjkMc0MeSwp7ynkE8yvQgLzVRMLmVZ8hH7usyLkgsuZdf9LLNofs+X
> siRZGG9DAgB4sGEREQSC+IYJBlGc/YpYs9VPSB5xS5vAGMEAjGAgD2uhGEK0IAEgSKvwqETUGAbh
> Z2dgiFUEH7lP+nE7+5D0qGKqj+8CKHR5R6DyL/BQ8h37i5Z6L8YODc8JFkgWtnydeOK4Z0Zi1C9Z
> y+y92+ITmq8VMr0r0l7HspEET6bxFzcMGMyqI5gxaoz4WWXOg5RM69rCUsJRPlClB5aKFORhRtiY
> E9Pq2NTrx1DJSggtqw2EbFr70FyaeF3RIodybV7FOiH4PAdwO/SE4oi7WECIkA8SIoFAWwKFL6q4
> uKofRQPFseZDu7CQ5KFUamAggO4bqa+S/ChrWwiOthE5VekPbQfghp0CIdXQhLE3c9l4+kEpHRI5
> AeZQ+VXiFLmCdBnOvoeR8QlKI9U51QoXah2g2VD50B8GnyKVvb/vQogQnygJ9gqsVEJpkBEZC5P5
> jUyQEYJpEu239ybQJ+CIwjGQEiRQP8PwpnCCCzP/V0/j682L4/lNxaP1qfj12/3b5bTR1BdBGEWC
> QsdhGYVBKqsapcB1sl0Ngo770ZUkBRUbSg6ppkH07myUGyerY/ClN+/9/6dRP0/Vm9FIwRn+T+5e
> y6pdIzbcfzosvGMPd1aJwwOHl//Ah8wQ//i7kinChIJcw/so
> _______________________________________________
> Mailing list: https://launchpad.net/~dolfin
> Post to : dolfin@xxxxxxxxxxxxxxxxxxx
> Unsubscribe : https://launchpad.net/~dolfin
> More help : https://help.launchpad.net/ListHelp
Attachment:
signature.asc
Description: Digital signature
Follow ups
References
