Calcola l'arricciatura di un campo vettoriale in Python e tracciarlo con matplotlib

Question

Ho bisogno di calcolare l'arricciatura di un campo vettoriale e tracciarlo con matplotlib. Un semplice esempio di ciò che sto cercando potrebbe essere:

Come posso calcolare e tracciare l'arricciatura del campo vettoriale nello quiver3d_demo.py nella galleria matplotlib?

Answer 1

È possibile utilizzare sympy.curl() per calcolare l'arricciatura di un campo vettoriale.

Esempio:

Supponiamo di avere:
F = (y z, -xy, z) = y z x - xy y + z z, quindi y sarebbe R[1], x è R[0] e z è R[2] mentre i vettori dei 3 assi sarebbero R.x, R.y, R.z e il codice per calcolare il campo vettoriale arricciatura è:

from sympy.physics.vector import ReferenceFrame 
from sympy.physics.vector import curl 
R = ReferenceFrame('R') 

F = R[1]**2 * R[2] * R.x - R[0]*R[1] * R.y + R[2]**2 * R.z 

G = curl(F, R) 

In tal caso G sarebbe pari a R_y**2*R.y + (-2*R_y*R_z - R_y)*R.z o, in altre parole,
G = (0, y , -2yz-y).

Per tracciarlo è necessario convertire il risultato sopra in 3 funzioni separate; u, v, w.

(esempio di seguito adattato dal matplotlib example on this link):

from mpl_toolkits.mplot3d import axes3d 
import matplotlib.pyplot as plt 
import numpy as np 

fig = plt.figure() 
ax = fig.gca(projection='3d') 

x, y, z = np.meshgrid(np.arange(-0.8, 1, 0.2), 
         np.arange(-0.8, 1, 0.2), 
         np.arange(-0.8, 1, 0.8)) 

u = 0 
v = y**2 
w = -2*y*z - y 

ax.quiver(x, y, z, u, v, w, length=0.1) 

plt.show() 

E il risultato finale è questo:

Answer 2

Per calcolare il ricciolo di una funzione vettoriale è anche possibile utilizzare numdifftools per differenziazione numerica automatica senza una deviazione attraverso la differenziazione simbolica. Numdifftools non fornisce una funzione curl(), ma calcola la matrice Jacobiana di una funzione di valore vettoriale di una o più variabili e fornisce le derivate di tutte le componenti di un campo vettoriale rispetto a tutte le variabili; questo è tutto ciò che è necessario per il calcolo del ricciolo.

import import scipy as sp 
import numdifftools as nd 

def h(x): 
    return sp.array([3*x[0]**2,4*x[1]*x[2]**3, 2*x[0]]) 

def curl(f,x): 
    jac = nd.Jacobian(f)(x) 
    return sp.array([jac[2,1]-jac[1,2],jac[0,2]-jac[2,0],jac[1,0]-jac[0,1]]) 

x = sp.array([1,2,3)] 
curl(h,x) 

Questo restituisce il valore del ricciolo a x: array([-216., -2., 0.]) Plotting è come suggerito sopra.