Mathematics made simple.

Followers

Showing posts with label Programming. Show all posts
Showing posts with label Programming. Show all posts

Generating Sexy Prime Pairs using Python


Sexy Primes, What are they?

They are prime numbers that differ by 6. For example, 5 & 11,  7 & 13. The name comes from the Latin word for six; sex. What we are going to do is, write a program to generate all the sexy prime pairs within a given interval of natural numbers. 

Sexy Prime Pairs

This is a program to generate all the sexy prime pairs below 10000.

Program

def check_prime(n):
    for i in range(2,n//2+1):
        if n%i == 0:
    return False
    return True

def primes_list(a,b):
    primes_list = []
    for i in range(a,b):
        if check_prime(i):
    primes_list.append(i)
    return (primes_list)

def sexy_list(primes_list):
    sexy_list = []
    for i in primes_list:
        for j in primes_list:
    if j-i == 6:
        sexy_list.append((i,j))
    return(sexy_list)

primes = primes_list(2,1000)
print(sexy_list(primes))

Output

[(5, 11), (7, 13), (11, 17), (13, 19), (17, 23), (23, 29), (31, 37), (37, 43), (41, 47), (47, 53), (53, 59), (61, 67), (67, 73), (73, 79), (83, 89), (97, 103), (101, 107), (103, 109), (107, 113), (131, 137), (151, 157), (157, 163), (167, 173), (173, 179), (191, 197), (193, 199), (223, 229), (227, 233), (233, 239), (251, 257), (257, 263), (263, 269), (271, 277), (277, 283), (307, 313), (311, 317), (331, 337), (347, 353), (353, 359), (367, 373), (373, 379), (383, 389), (433, 439), (443, 449), (457, 463), (461, 467), (503, 509), (541, 547), (557, 563), (563, 569), (571, 577), (587, 593), (593, 599), (601, 607), (607, 613), (613, 619), (641, 647), (647, 653), (653, 659), (677, 683), (727, 733), (733, 739), (751, 757), (821, 827), (823, 829), (853, 859), (857, 863), (877, 883), (881, 887), (941, 947), (947, 953), (971, 977), (977, 983), (991, 997)]

Monte Carlo Integration using Python


 Monte Carlo Integration

Monte Carlo integration is a  powerful method of computing the value of complex integrals using probabilistic techniques. This technique uses random numbers to compute the definite integral of a function. Here we are going to use the python programs written in the previous post to generate pseudorandom numbers and approximate the value of the definite integral of a function.
Consider a function to be integrated, as shown below:

That is we need to evaluate the definite integral $\int_{a}^{b} f(x)\mathrm{d}x$. The integral is just the area under the curve. The width of the interval $(b- a)$ times the average value of the function is also the value of the integral, that is,

$\int_{a}^{b}f(x)\mathrm{d}x = (b - a)f_{average} = (b-a)\langle f \rangle$

So if we had some independent way of calculating the average value of the integrand, then we could evaluate the integral. This is where random numbers come in.
Imagine that we had a list of random numbers, $x_j$ uniformly distributed between $a$ and $b$. To calculate the function average, we simply evaluate $f(x_j)$ at each of the randomly selected points, add all of them together and divide the sum by the number of points.

$$\langle f \rangle _N = \frac{1}{N}\sum_{i=1}^{N}f(x_i)$$

As the number of points used in calculating the average increases, $\langle f \rangle_N$ approaches the true value, $\langle f \rangle$. Therefore, as a numerical approximation, we can write

$$\int_{a}^{b}f(x)\mathrm{d}x = \frac{(b-a)}{N}\sum_{i=1}^{N}f(x_i)$$

Alternatively, we can look at the Monte Carlo integration method in the following way:
  • Find some value $M$ such that $f(x) < M$ over the integral $[a,b]$.
  • Select random numbers $x$ and $y$ from the intervals $[a,b]$ and $[0,M]$ respectively.
  • Determine if $y > f(x)$ or $y \leq f(x)$.
  •  Repeat the process $N$ times, keeping track of the number of times $y \leq f(x)$ or under the curve (successes); call the total number of successes $S$.
The estimated probability of success is

$$\frac{S}{N} = \frac{Area \ under \ curve}{Total \ area \ inside \ rectangle} = \frac{\int_{a}^{b}f{x}\mathrm{d}x}{M(b-1)}$$

Thus the value of the integral can be calculated from the formula

$$\int_{a}^{b}f(x)\mathrm{d}x = M(b-a)\frac{S}{N}$$

The python program that implements the Monte Carlo method of integration is given below. Here the random numbers used in the integration is generated using Lehmer's method mention in the previous post.

Program

from sympy import *
def pseudo(a,m,x0):
    list1 = []
    for i in range(100):
        x0 = (a*x0)%m
        list1.append(x0)
    return(list1)

def integrate(f,var):
    f = sympify(f)
    var = Symbol(var)
    rand_values1 = pseudo(61,101,18)
    rand_values2 = pseudo(63,107,18)
    f_values1 = []
    f_values2 = []
    x = []
    y = []
    for i in rand_values1:
        x.append(a + (((b-1)/101)*i))
    i = a
    #to find maximum M
    while i <= b:
        f_values1.append(f.subs({var:i}))
        i += 0.01
    #to check whether y lies below or above the curve
    for i in x:
        f_values2.append(f.subs({var:i}))
    M = int(max(f_values1)) + 1
    for i in rand_values2:
        y.append((M/107)*(i))
    count = 0
    for i in range(len(y)):
        if y[i] <= f_values2[i]:
            count += 1
    return (b-a)*M*(count/100)

f = input("Enter the function: ")
var = input("Enter the variable: ")
a = int(input("Enter the lower limit: "))
b = int(input("Enter the upper limit: "))
print("The value of the definite integral is", integrate(f,var))


The value of some definite integrals obtained using the above program and their actual values are given below in the table:


Generating Pseudorandom Numbers Using Congruence Theory In Python


Pseudorandom Numbers

A sequence of random numbers, chosen from the interval between 1 and some fixed integer $M$, is a sequence of numbers $x_0, x_1, \ldots $ such that for each $i$, the values of $x_i$  has a $1/M$ chance of being any given number between 1 and $M$, independent of what values $x_0, x_1, \ldots ,x_{i-1}$ took on.
For example, $M$ could be 6 and $x_i$ could be the number of dots showing on the $i^{th}$ throw of a fair die.
Random numbers are useful in many contexts related to computing. For example, cryptography requires numbers that attackers can't guess. We can't just use the same numbers over and over. Thus it is crucial to generate these numbers in a very unpredictable way so that attackers can't guess them.
We generally group the random numbers generated by computers into two types, depending on how they're generated: True random numbers and pseudorandom numbers. To generate a true random number, the computer measures some type of physical phenomenon that takes place outside the computer system. For example, the computer could measure the radioactive decay of an atom, the atmospheric noise or simply use the exact time one press keys on a keyboard as a source of unpredictable data.
Since a computer is not a random device, but rather a deterministic machine, it cannot compute random numbers on its own. It needs some kind of data from the physical world for that as discussed above. To circumvent this problem, a strategy adopted by computer scientists has been for the computer to generate sequences of numbers that look random, even though in fact they are not random. Such numbers are called pseudorandom numbers.
Among all the methods used to compute pseudorandom numbers, the most widely used and understood is the multiplicative congruential method developed by D.H. Lehmer in 1948.
The method developed by Lehmer is as follows:
  • Pick a number $M$, the modulus.
  • Pick a starting number $X_0$, the seed value.
  • Pick a multiplier $A$.
  • Generate the sequence as follows:
    $X_1$ is the remainder obtained on dividing $AX_0$ by $M$; that is $AX_0 \pmod M$
    $X_2$ is the remainder obtained on dividing $AX_1$ by $M$; that is $AX_1 \pmod M$
    $\vdots$
    $X_{n+1}$ is the remainder obtained on dividing $AX_n$ by $M$; that is $AX_n \pmod M$
Thus each of $X_0, X_1, X_2, \ldots$ is a number between 0 and $M$, and each $X_n$ satisfies
$X_{n+1} \equiv AX_n \pmod M$
To implement Lehmer's method in practice to get a good sequence of pseudorandom numbers, one should use a large modulus $M$. A large modulus $M$ commonly used is $M = 2^{31} - 1$ and the multiplier $A$ widely used is $A =7^5$. The python code that implement the Lehmer's method of pseudorandom number generation is given below for $M = 2^{31} - 1, A = 7^{5}$ and $X_0 = 621317$. The program displays the first 100 numbers in the pseudorandom sequence.  

Program

def pseudo(a,m,x0):
    list1 = []
    for i in range(100):
        x0 = (a*x0)%m
        list1.append(x0)
    print(list1)
m,a,x0 = 2**31-1,7**5,621317
pseudo(a,m,x0)

Output

[1852540231, 1425748211, 927649051, 266322937, 733681811, 139096403, 1331037285, 406498196, 869699065, 1258483973, 773694908, 476836171, 1924039040, 515388754, 1337240127, 1578448634, 1120700247, 29983492, 1421376646, 469200094, 286028074, 1205437732, 431235926, 24899657, 1876707681, 1733671078, 751685450, 2078546496, 1014472523, 1367020528, 1733958490, 1287251640, 1088053602, 1093634609, 404338790, 1083784422, 208486700, 1490138643, 805881587, 272471080, 986306156, 421292699, 412807934, 1690766928, 1216141792, 2093229645, 833538361, 1243403946, 726751465, 1772371766, 506603625, 1861948667, 639542185, 629850060, 943062357, 1619719239, 1118540501, 238354469, 966558828, 1387916288, 741678702, 1398857326, 2091594373, 1266809268, 1110490918, 240482749, 229338789, 1911364005, 86956562, 1190057574, 1782441707, 100893899, 1359163010, 669155931, 131872978, 186017542, 1808122009, 65516566, 1625297498, 383059046, 2064896063, 1372395321, 1873791267, 2109624861, 1510026857, 59645353, 1732068369, 1732242698, 367222907, 47396471, 2023538707, 2064014657, 1590990208, 1453537059, 1970865988, 1556888988, 1692466268, 1859661761, 858218689, 1581332771]

Now we can visualize the first 1000 pseudorandom numbers generated using the above algorithm by a python library called matplotlib.

Program

import matplotlib.pyplot as plt
def pseudo(a,m,x0):
    x = []
    y = []
    for i in range(1000):
        x0 = (a*x0)% m
        x.append(i)
        y.append(x0)
    return(x,y)
m,a,x0 = 2**31 - 1, 7**5, 621317
x,y = pseudo(a,m,x0)
plt.plot(x,y,"o")
plt.show()

Output

The Scatter Plot of the First 1000 Pseudorandom Number

Note. For the sake of programming, $X_0$ is removed from the first position of the list. But we can see that $X_0$ appears again later at some point in the list. In fact,  after $X_0$ each element in the list starts repeating itself. 

Making an Encryption Application in Python Using the RSA Algorithm


 The RSA Cryptosystem

RSA algorithm is an asymmetric cryptography algorithm. Asymmetric actually means that it works on two different keys i.e. Public Key and Private Key. As the name describes that the public key is given to everyone and the private key is kept private.
The Idea: The idea of the RSA algorithm is based on the fact that it is difficult to factorize a large integer. The public key consists of two numbers where one number is the product of two large prime numbers. And private key is also derived from the same two prime numbers. So if somebody can factorize the large number, the private key is compromised. Therefore encryption strength totally lies on the key size and if we double or triple the key size, the strength of encryption increases exponentially. RSA keys can be typically 1024 or 2048 bits long, but experts believe that 1024 bit keys could be broken in the near future. But till now it seems to be an infeasible task.

The mechanism behind the RSA algorithm




Now we have our public (n = 3127, e = 3) and private (d = 2011) keys. 

Below is a Python implementation of the RSA algorithm:

Here a python package called tkinter is used to make the graphical user interface (GUI) for the application we are developing. 

Program

import tkinter as tk
from tkinter import ttk
letter = [ ”a” , ”b” , ”c” , ”d” , ”e” , ”f” , ”g” , ”h” , ”i” , ”j” , ”k” , ”l” , ”m” , ”n” , ”o” , ”p” , ”q” , ”r” , ”s” , ”t” , ”u” , ”v” , ”w” , ”x” , ”y” , ”z” , ” ” , ”A” , ”B ” , ”C” , ”D” , ”E” , ”F” , ”G” , ”H” , ”I” , ”J” , ”K” , ”L” , ”M” , ”N” , ”O” , ”P” , ”Q” , ”R” , ”S” , ”T” , ”U” , ”V” , ”W” , ”X” , ”Y” , ”Z” ] 
n = 53*59
phi = 52*58
e = 3
d = 2011
a = []
b = []
c = []
f = []
def encryption():
    msg = entry1.get()
    for i in msg:
        pos = letters.index(i)
        a.append(pos)
    for j in a:
        b.append((j**e)%n)
    label2.configure(text = ”The Encrypted text”)
    label3.configure(text = b)
def decryption():
    msg = entry1.get()
    g = " "
    l = list(map(int,msg.split(" ")))
    for k in l:
        c.append((k**d)%n)
    for r in c:
        f.append(letters[r])
    for j in f:
        g += j
    label2.configure(text = "The Decrypted text")
    label3.configure(text = g)
root = tk.Tk()
root.title("RSA Algorithm")
root.geometry("500x300")
label1 = ttk.Label(text = "Enter the Message to be Encrypted or Decrypted")
label1.place(x =130, y = 50)
label2 = ttk.Label(text = " ")
label2.place(x = 100,y = 180)
label3 = ttk.Label(text = " ")
label3.place(x = 100, y =200)
entry1 = ttk.Entry()
entry1.place(x=175,y=90)
button1 = ttk.Button(text = "Encrypt",cursor="hand2",command=encryption)
button1.place(x=100,y=140)
button2 = ttk.Button(text = "Decrypt" , cursor = "hand3" , command = lambda:entry1.delete(0,tk.END))
button3.place(x=200,y=140)
root.mainloop()

Output

Encryption

Decryption

Modular Designs With Python (Congruence Theory)


Modular Designs

Congruence theory can be used to generate beautiful designs. It is really a fun way to study congruences. Two such designs are created using python. $m$-pointed stars and $(m,n)$ residue designs.

$m$ - Pointed Stars

To construct an $m$-pointed star, mark $m$ equally spaced points on a large circle, and label them with the least residue 0 through $(m-1)$ modulo $m$. Choose a least residue $i$ modulo $m$, where ($i,m$) = 1. Join each point $x$ with the point  $x + r$ modulo $m$. Now colour in the various regions inside the circle with some solid colours. We should get a nice $m$-pointed star.
A python program can be written to draw an $m$-pointed star, where the value of m and $r$ can be chosen by the user.

Program

import turtle
from math import *
m = int(input("Enter the value of m: "))
r = int(input("Enter the value of r: "))
t = turtle.Turtle()
t.speed(3)
t.penup()
t.color("black","cyan")
t.setposition(0,-100)
t.pendown()
a=[]
for i in range(m):
    t.circle(100,360/m)
    t.dot(5,"blue")
    s = t.position()
    t.write(i)
    a.append(s)
n = len(a)
t.penup()
t.setposition(a[0])
t.pendown()
i = 0
t.begin_fill()
while i <= m*r:
    t.setposition(a[i%m])
    i+=r
t.end_fill()

Output

A seven-pointed (r = 3) star 
 
A thirteen-pointed (r = 5) star



$(m,n)$ Residue Designs

To construct an $(m,n)$ residue design, where $1 \leq n < m$ and $(m,n) = 1$, select $m-1$ equally spaced points on a large circle. Label them 1 through $m -1$, and join each $x$ to point $nx$ modulo $m$. Then color in the various regions formed in a systematic way to create exciting designs.
The python program to do the same is given below. The program asks the user to input the values of $m$ and $n$.

Program

import turtle
m = int(input("Enter the value of m: "))
n = int(input("Enter the value of n: "))
t = turtle.Turtle()
t.speed(3)
t.color("black","cyan")
a = []
for i in range(1,m):
    t.circle(100,360/(m-1))
    t.dot(5,"blue")
    t.write(i)
    a.append(t.position())
t.penup()
t.setposition(a[0])
t.pendown()
a.append((0,0))
i = 0
t.begin_fill()
while i < m-1:
    t.goto(a[(((i+1)*n)%m)-1])
    t.penup()
    i += 1
    t.goto(a[i])
    t.pendown()
t.end_fill()

Output

The (19,9) reside design
The (21,10) residue design

Visualization of Collatz Conjecture using Python


The Collatz Conjecture



What do you think of this image? Looks like a feather. Right? Well, this feather is drawn by using a well-defined mathematical sequence. To understand how this feather is drawn we need to understand a 60-year-old problem in mathematics.

The Problem

The Collatz conjecture (a.k.a the hailstone problem or the $3n + 1$ problem) was proposed by Lother Collatz in 1937. Although the problem on which the conjecture is based is really simple that even a fourth-grader can easily understand it, the behaviour of the conjecture makes it exceedingly difficult to prove(or disprove). 
First of all, lets define the Collatz map. Let $T : \mathbb{Z}^+ \rightarrow \mathbb{Z}^+$ be defined by

The conjecture is that for every positive integer $n$, there exists a $k$ such that  $T^k(n) = 1$.

Let us consider an example.

As mentioned above, the conjecture states that this is true for any positive integer n. But what makes this problem interesting is that even after 60  years since its proposal nobody has been able to prove (or disprove) it. Mathematicians say the "Mathematics is yet not ready to tackle such problems." So far the conjecture has been checked for all starting values up to $87 \times 2^{60}$. All initial values tested so far eventually end up in 1.

Program

#Program to check the Collatz conjecture for a given range of numbers.

def collatz(x):
    a = [x]
    while x!=1:
        if x % 2 == 0:
            x = x/2
        else:
            x = 3*x + 1
        a.append(x)
    if a[len(a)-1] == 1:
        return("Yes")
    else:
        return("No")
for i in range(1,100000):
    print(i,collatz(i))

This program checks the conjecture for the first 100000 positive integers, all of which returns to be true.

Visualization

Now, what we are really interested is on the visualization of the problem. We try to visualize the problem with the aid of Python programs and try to find some patterns that may arise in the figures.
  • Stopping time
    Python program to compute the stopping time $k$ for each $n$ and plot them against $n$ in a graph.

    Program
    #Program to compute the stopping time for each n and to plot it in a graph.
    import matplotlib.pyplot as plt
    def collatz(x):
        a = [x]
        while x != 1:
            if x % 2 == 0:
                x = x/2
            else:
                x = 3*x + 1
            a.append(x)
        return(a)
    x=[]
    y=[]
    for i in range(1,10000):
        a = collatz(i)
        x.append(i)
        y.append(len(a))
    plt.plot(x,y,"o")
    plt.show()

    We can see a pattern in the obtained graph.


  • The highest value in the iteration
    For each $n$, the composition of the function $T$ yields an upper bound. For example, 16 in the case of n = 6. Let's plot the graph of this largest value obtained in the iteration for each n and try to figure out some patterns from the graph.|
    Program
    #Highest value in the iteration
    import matplotlib.pyplot as plt
    def collatz(x):
        a = [x]
        while x != 1:
            if x % 2 == 0:
                x =x/2
            else:
                x = 3*x + 1
            a.append(x)
        return(a)
    b = []
    c = []
    for i in range(1,1000):
        c.append(max(collatz(i)))
        b.append(i)
    plt.plot(b,c,"o")
    plt.show()
From the graph, it is evident that for the first 1000 positive integers the highest value obtained in the iteration for almost all $n$'s is bounded in the y-axis. For some $n$ they wander off to large numbers.

The Feather

Now let's get back to the feather we discussed earlier. Let's make a directed graph for the conjecture. It somewhat looks like this:


Now starting from 1 move upward along the graph. An edge incident to an even number should be curved slightly towards the right and the one incident to an odd number should be curved slightly towards the left. Lose the numbers on the graph and colour the graph according to your imagination. This yields the feather.

Ultimate Theorem. Powered by Blogger.

About Me

My photo
Mathematics has always fascinated me. I love the subject since childhood. This love towards the field helped me in completing a Masters Degree in Mathematics. As much as I love the subject, I love teaching it too.

Local Linear Approximation

You might have seen the famous approximation $\sin(x) \approx x$ for $x$ near 0. We can use derivatives to approximate non-linear functions ...

Search This Blog

ultimate theorem