Plotting 3D level sets in Paraview

December 30, 2017 Leave a comment

Surfaces can be represented as certain level-sets for some 3D functions. Given a set of points with values attached to them a level set associated to a certain number will separate the points into two sets: with values higher or lower than the number considered. It is nice to have good looking plots when working with the level-set method in shape optimization. Paraview is a nice, open source framework, which has the right tools in order to produce high quality plots.

I’ll briefly describe below how to use Paraview to make some nice pictures of level-sets. First of all, you’ll need your data in some format that Paraview can understand. I use vtk file format for which there is a nice automated interface in the software I use for the optimization (FreeFem++). In the vtk file you need to have a set of points and a scalar value attached to them.

If you want to create level-sets associated to certain values, follow the steps below:

    1. Load your file (containing points with at least a scalar value) and click Apply.
    2. Next, go to Filters/Alphabetical and select “Cell Data to Point Data” (if you forget to do this, you’ll get a rough surface where you see the discretization). Click Apply.
    3. Then apply the Contour filter (by clicking the button or going in Filters/Alphabetical). You’ll have to select the field for which the contours will be drawn and then put in the values of the level-sets you want to see. Click Apply. An example can be seen below.lvlview

If your level-set cuts the boundary of the domain, Paraview will draw a hole there. If you want to have a closed region instead, you need to use the IsoVolume filter instead of the Contour one. The difference is that you need to specify two values and Paraview will draw the surface enclosing the points corresponding to these values. Many other features are directly available: you can color the level set following another scalar value, you can set the lighting, etc. You can also symmetrize your geometry using the Reflect filter. Below you can see a result built from my work.


You can also create animations in a pretty  straightforward way. Just go to View and select the Animation box. Then you’ll see the animation options. Add a Camera object with Orbit field selected. You’ll be presented with multiple options, like the center of rotation, direction of the vertical and initial position. Once everything is set, click the play button to see the animation. Then go to File/Save Animation to save it to a file.


I heard that Paraview could to many things when dealing with visualization aspects, but I hesitated to use it until now since the interface is not straightforward. The use of Filters is not clear in the beginning, but after playing with some examples, everything becomes really easy to use. The next step is to automatize all this using scripts.

Happy New Year!


Putnam 2017 A3 – Solution

December 4, 2017 Leave a comment

Problem A3. Denote {\phi = f-g}. Then {\phi} is continuous and {\int_a^b \phi = 0}. We can see that

\displaystyle I_{n+1}-I_n = \int_a^b (f/g)^n \phi = \int_{\phi\geq 0} (f/g)^n \phi+ \int_{\phi<0} (f/g)^n \phi

Now note that on {\{ \phi>=0\}} we have {f/g \geq 1} so {(f/g)^n \phi \geq \phi}. Furthermore, on {\{\phi<0\}} we have {(f/g)^n <1} so multiplying with {\phi<0} we get {(f/g)^n \phi \geq \phi}. Therefore

\displaystyle I_{n+1}-I_n \geq \int_{\phi \geq 0} \phi + \int_{\phi<0} \phi = \int \phi = 0.

To prove that {I_n} goes to {+\infty} we can still work with {I_{n+1}-I_n}. Note that the negative part cannot get too big:

\displaystyle \left|\int_{ \phi <0 } (f/g)^n \phi \right| \leq \int_{\phi<0} |\phi| \leq \int_a^b |f-g|.

As for the positive part, taking {0<\varepsilon< \max_{[a,b]} \phi} we have

\displaystyle \int_{\phi\geq 0} (f/g)^n \phi \geq \int_{\phi>\varepsilon}(f/g)^n \varepsilon.

Next, note that on {\{ \phi \geq \varepsilon\}}

\displaystyle \frac{f}{g} = \frac{g+\phi}{g} \geq \frac{g+ \varepsilon}{g}.

We would need that the last term be larger than {1+\delta}. This is equivalent to {g\delta <\varepsilon}. Since {g} is continuous on {[a,b]}, it is bounded above, so some delta small enough exists in order for this to work.

Grouping all of the above we get that

\displaystyle I_{n+1}-I_n \geq \int_{\phi \geq 0} (f/g)^n \phi \geq \int_{\phi>\varepsilon} \varepsilon (1+\delta)^n.

Since {|\phi>\varepsilon|>0} we get that {I_{n+1}-I_n} goes to {+\infty}.

Putnam 2017 A2 – Solution

December 4, 2017 Leave a comment

Problem A2. We have the following recurrence relation

\displaystyle Q_n = \frac{Q_{n-1}^2-1}{Q_{n-2}},

for {n \geq 2}, given {Q_0=1} and {Q_1=x}. In order to prove that {Q_n} is always a polynomial with integer coefficients we should prove that {Q_{n-2}} divides {Q_{n-1}^2-1} somehow. Recurrence does not seem to work very well. Also, root based arguments might work, but you need to take good care in the computation.

A simpler idea, which is classic in this context, is to try and linearize the recurrence relation. In order to do this, let’s write two consecutive recurrence relations

\displaystyle Q_nQ_{n-2} +1 = Q_{n-1}^2

\displaystyle Q_n^2 = Q_{n+1}Q_{n-1}+1

We add them and we obtain the following relation

\displaystyle \frac{Q_n}{Q_{n-1}} = \frac{Q_{n+1}+Q_{n-1}}{Q_n+Q_{n-2}},

which leads straightforward to a telescoping argument. Finally, we are left with a simple linear recurrence with integer coefficient polynomials, and the result follows immediately.

Putnam 2017 – Problem A1

December 4, 2017 Leave a comment

Problem A1. We have {n^2 \in S \Rightarrow N \in S} and {n \in S \Rightarrow (n+5)^2 \in S}. Therefore {n \in S \Rightarrow n+5 \in S}.

Next, let’s note what elements cannot be in {S}. Note that taking square roots and squaring {n+5} cannot change a non-zero remainder modulo {5} into a zero remainder. Therefore, starting from {2} one could never get a multiple of {5} following the allowed operations. Thus we can safely say that multiples of {5} are not in the minimal set {S}.

Furthermore, {1} could only be obtained as a square root of itself with the allowed operations. Starting from {2}, one could never get below by performing square roots or {n \mapsto (n+5)^2}. Therefore, the minimal set {S} does not contain {1} and multiples of {5}.

Now, we show that it contains all the rest. The general idea is as follows: it is enough to find which is the smallest element in a class of remainders modulo {5} to deduce that all larger elements are there (recall the operation {n \in S \Rightarrow n+5\in S}). Now in order to obtain small elements of {S}, one would need to take successive square roots. So if we prove that for some {a} we have {a^{2^n}\in S} for some {n} then we get that {a \in S}.

Now let’s start from the beginning. We have {2 \in S} so {49 = (2+5)^2 \in S}. Since {49+5k} is in {S} for every {k}, we get that all squares of the form {5k+4} greater than {49} are in {S}. Moreover, {(49+5)^2 = 2916} so all numbers of the form {5k+1} greater than {2916} are in {S}. Since {81^2 =6561} it follows that {6561 \in S \Rightarrow 81 \in S \Rightarrow 9 \in S \Rightarrow 3 \in S}. Moreover, {6^{16}} ends in {6} and is greater than {2916} so {6 \in S}. Next, we have {4^8 = 16^4} which ends in {6} and is greater than {2916} so it is also in {S}. Therefore {4 \in S}.

Finally, we have that {n \in S \Rightarrow n+5 \in S}, {\{2,3,4,6\} \subset S} and {1 \notin S}, {5k \notin S}. This means that the minimal set {S} is {\Bbb{Z}_+^* \setminus\{1\} \setminus \{5k: k \in \Bbb{Z}_+\}}.

Categories: Uncategorized

Putnam 2017 – Problems

December 4, 2017 Leave a comment

Source: Art of Problem Solving forum

Problem A1. Let {S} be the smallest set of positive integers such that

  • a) {2} is in {S,}
  • b) {n} is in {S} whenever {n^2} is in {S,} and
  • c) {(n+5)^2} is in {S} whenever {n} is in {S.}

Which positive integers are not in {S?}

(The set {S} is “smallest” in the sense that {S} is contained in any other such set.)

Problem A2. Let {Q_0(x)=1}, {Q_1(x)=x,} and

\displaystyle Q_n(x)=\frac{(Q_{n-1}(x))^2-1}{Q_{n-2}(x)}

for all {n\ge 2.} Show that, whenever {n} is a positive integer, {Q_n(x)} is equal to a polynomial with integer coefficients.

Problem A3. Let {a} and {b} be real numbers with {a<b,} and let {f} and {g} be continuous functions from {[a,b]} to {(0,\infty)} such that {\int_a^b f(x)\,dx=\int_a^b g(x)\,dx} but {f\ne g.} For every positive integer {n,} define

\displaystyle I_n=\int_a^b\frac{(f(x))^{n+1}}{(g(x))^n}\,dx.

Show that {I_1,I_2,I_3,\dots} is an increasing sequence with {\displaystyle\lim_{n\rightarrow\infty}I_n=\infty.}

Problem A4. A class with {2N} students took a quiz, on which the possible scores were {0,1,\dots,10.} Each of these scores occurred at least once, and the average score was exactly {7.4.} Show that the class can be divided into two groups of {N} students in such a way that the average score for each group was exactly {7.4.}

Problem A5. Each of the integers from {1} to {n} is written on a separate card, and then the cards are combined into a deck and shuffled. Three players, {A,B,} and {C,} take turns in the order {A,B,C,A,\dots} choosing one card at random from the deck. (Each card in the deck is equally likely to be chosen.) After a card is chosen, that card and all higher-numbered cards are removed from the deck, and the remaining cards are reshuffled before the next turn. Play continues until one of the three players wins the game by drawing the card numbered {1.}

Show that for each of the three players, there are arbitrarily large values of {n} for which that player has the highest probability among the three players of winning the game.

Problem A6. The {30} edges of a regular icosahedron are distinguished by labeling them {1,2,\dots,30.} How many different ways are there to paint each edge red, white, or blue such that each of the 20 triangular faces of the icosahedron has two edges of the same color and a third edge of a different color?

Problem B1. Let {L_1} and {L_2} be distinct lines in the plane. Prove that {L_1} and {L_2} intersect if and only if, for every real number {\lambda\ne 0} and every point {P} not on {L_1} or {L_2,} there exist points {A_1} on {L_1} and {A_2} on {L_2} such that {\overrightarrow{PA_2}=\lambda\overrightarrow{PA_1}.}

Problem B2. Suppose that a positive integer {N} can be expressed as the sum of {k} consecutive positive integers

\displaystyle N=a+(a+1)+(a+2)+\cdots+(a+k-1)

for {k=2017} but for no other values of {k>1.} Considering all positive integers {N} with this property, what is the smallest positive integer {a} that occurs in any of these expressions?

Problem B3. Suppose that

\displaystyle f(x) = \sum_{i=0}^\infty c_ix^i

is a power series for which each coefficient {c_i} is {0} or {1}. Show that if {f(2/3) = 3/2}, then {f(1/2)} must be irrational.

Problem B4. Evaluate the sum

\displaystyle \sum_{k=0}^{\infty}\left(3\cdot\frac{\ln(4k+2)}{4k+2}-\frac{\ln(4k+3)}{4k+3}-\frac{\ln(4k+4)}{4k+4}-\frac{\ln(4k+5)}{4k+5}\right)

\displaystyle =3\cdot\frac{\ln 2}2-\frac{\ln 3}3-\frac{\ln 4}4-\frac{\ln 5}5+3\cdot\frac{\ln 6}6-\frac{\ln 7}7-\frac{\ln 8}8-\frac{\ln 9}9+3\cdot\frac{\ln 10}{10}-\cdots.

(As usual, {\ln x} denotes the natural logarithm of {x.})

Problem B5. A line in the plane of a triangle {T} is called an equalizer if it divides {T} into two regions having equal area and equal perimeter. Find positive integers {a>b>c,} with {a} as small as possible, such that there exists a triangle with side lengths {a,b,c} that has exactly two distinct equalizers.

Problem B6. Find the number of ordered {64}-tuples {\{x_0,x_1,\dots,x_{63}\}} such that {x_0,x_1,\dots,x_{63}} are distinct elements of {\{1,2,\dots,2017\}} and

\displaystyle x_0+x_1+2x_2+3x_3+\cdots+63x_{63}

is divisible by {2017.}

Categories: Uncategorized

A hint for Project Euler Pb 613

November 18, 2017 Leave a comment

The text for Problem 613 can be found here. The hint is the following picture 🙂


Metapost – or how to code an image

November 9, 2017 Leave a comment

What software do you use when you need to draw a nice image? I often need to draw things related to research and I never managed to efficiently use a software which uses the mouse or touchpad to draw and modify things. Moreover, if you need to add some mathematical text to the figure things get even more complicated. For me it is more natural to use code to generate graphics, if this is possible. Using something like Matlab to draw is possible. The advantage is that once you have a working code which produces what you want, you can easily modify it. If you have an image where you need to repeat things, using loops can facilitate the job (programmers will understand…).

This is were Metapost comes into play. I found this a long time ago while searching for a tool to build nice graphics for my master thesis. Being used to LaTeX for typesetting math, it was a natural way to draw using code. I do not claim it is the best or the simplest way, but for me it works. Most importantly, it allows me to build high quality, vectorized graphics, which are memory efficient. The advantage of vector graphics is that you can zoom as much as you want and you’ll never see the pixels.

There are a bunch of places where you can learn Metapost. I’ll put here some references I use a lot. The easiest way to start drawing in Metapost is to take a look at some existing codes and modify them. You can find lots of examples in this tutorial. If you have a linux system, you can compile metapost .mp files to obtain high quality pdfs using the command mptopdf. If not, then you can use the online Metapost editor found here. In any case, here are some of my codes for some drawings in Metapost.

Here is the code for my website logo: logo If you want to see the lossless vectorized pdf click to see the following file: logo-0    You can see that it has a border which changes from one color to another. This was done using a loop and a parameter to vary the color as we go along the boundary.

% copy from here to use the online previewer
u:=25; % 25 = 25bp = 25 PostScript points = 25/72 in
wi:=10; % width in units u
he:=7; % height in units u
hoehe:=he*u; % height
%for i=0 upto he:
%draw (0, i*u)--(breite, i*u) withcolor .7white;
%for j=0 upto wi:
%draw (j*u, 0)--(j*u, hoehe) withcolor .7white;
path p,q;
pickup pensquare scaled 20;
draw p;
numeric c,d,detail;
color a,b,co;
for i=1 upto (detail/2):
q := subpath (h*(i-1)/detail,h*i/detail) of p;
draw q withcolor (2*i/detail)[a,b];
for i=detail downto (1+detail/2):
q := subpath (h*(i-1)/detail,h*i/detail) of p;
draw q withcolor (2*(detail-i)/detail)[a,b];
label ("B",(1/15)*(2.6u,2.7u)) scaled 15 withcolor co ;
label ("2" infont defaultfont scaled 8,(4.9u,4.1u))
withcolor co;
% end copy here for online previewer

Here’s another example of a figure where I needed a grid of disks aligned on top of another figure. Using loops in Metapost allowed me to get what I wanted: cm

For the high quality pdf click here: cioranescu-murat-0  and see the code below

u:=25; % 25 = 25bp = 25 PostScript points = 25/72 in
wi:=10; % width in units u
he:=7; % height in units u
hoehe:=he*u; % height

path p,pa;
draw p;
fill p withcolor .8white;
pair a,b;
numeric c,d;
a:= point c of p;
b:= point d of p;
pickup pencircle scaled 1.5;
%draw a;
%draw b;
path q,r,s;
q:= subpath (c,d) of p;
%draw q withcolor red;
%r:= b..((a+b)/2 shifted (.5(a-b) rotated -90))..a;
%draw r;
%fill buildcycle(r,q) withcolor .5red;
%fill pa withcolor .9999blue;
%label.rt(btex $\Omega$ etex,.5(u,3u))  scaled 2;

draw (-u,9u)-- (10u,9u)--(10u,-2u)--(-u,-2u)--cycle;
for i=0 upto 9:
for j=-1 upto 8:
    draw fullcircle scaled 0.3u shifted (i*u, j*u);
    fill fullcircle scaled 0.3u shifted (i*u, j*u) withcolor white;

There are many ways today to draw what you want. Metapost is one of the tools available, if you like coding.

%d bloggers like this: