In [2]:
import numpy as np
import plotly
import matplotlib.pyplot as plt

plotly.offline.init_notebook_mode(connected=True)  # pour export html

def point(x, y, z):
    return np.array([x, y, z])

I Géométrie

In [3]:
# Q1
def vec(A, B):
    return np.array(B) - np.array(A)
In [4]:
# Q2 
def ps(v1, v2):
    return np.sum(v1 * v2)
In [5]:
# Q3
def norme(v):
    return np.sqrt(np.sum(v**2))
In [6]:
# Q4
def unitaire(v):
    return v / norme(v)
In [7]:
# Q5
def pt(r, t):
    S, u = r
    return S + t * u

def dir(A, B):
    return unitaire(vec(A, B))

def ra(A, B):
    return A, dir(A, B)
In [8]:
# Q6
def sp(A, B):
    return (A, norme(B - A))
In [9]:
# Q8
def intersection(r, s):
    S, u = r
    C, r_s = s
    CS = vec(S, C)
    b = 2*ps(u, CS)
    delta = b**2 - 4*(norme(CS)**2 - r_s**2)
    if delta < 0:
        return None
    t1, t2 = -b - delta**0.5, -b + delta**0.5
    t = min(abs(t1), abs(t2))/2
    return S + t*u, t

II Optique

Q 9. P est le premier point d'intersection d'un rayon partant de S et allant en P. (ou : Produit scalaire positif avec la normale au plan tangent)

In [10]:
# Q10
def au_dessus(s, P, src):
    return ps(vec(P, src), vec(s[0], P)) > 0
In [11]:
# Q11
def visible(obj, j, P, src):
    for i, o in enumerate(obj):
        inter = intersection(ra(src, P), o)
        if i != j and inter != None and inter[1] < norme(vec(src, P)):
            return False
    return au_dessus(obj[j], P, src)
In [12]:
# Q12
def couleur_diffusée(r, Cs, N, kd):
    cos = ps(-r[1], N)
    return cos*(Cs*kd)
In [13]:
# Q13
def rayon_réfléchi(s, P, src):
    C, r = s
    u = dir(src, P)
    N = dir(C, P)
    return P, 2*N + u

Digression : affichage d'un rayon avec Plotly

In [14]:
import plotly.graph_objects as go

def line(P, Q):
    return go.Scatter3d(x=[P[0], Q[0]], y=[P[1], Q[1]], z=[P[2], Q[2]])

theta, phi = np.mgrid[0:2*np.pi:100j, 0:2*np.pi:100j]
s = (point(0, 0, 0), 1)
src = point(3, 3, 2)
P = point(1, 0, 0)
_, u = rayon_réfléchi(s, P, src)

fig = go.Figure(data=[ 
    line(src, P),
    go.Surface(x=np.cos(theta)*np.cos(phi), y=np.sin(theta)*np.cos(phi), z=np.sin(phi)),
    line(P, P + 3*u),
])
fig.show(renderer="notebook")