Commit 62945ae1 authored by Jeremy BLEYER's avatar Jeremy BLEYER

Add simple file

parent 513729ad
#!/usr/bin/env python3
# -*- coding: utf-8 -*-
import numpy as np
from numba import jit, float64, int64, prange, generated_jit
import matplotlib.pyplot as plt
......@@ -22,12 +24,11 @@ def rec_fill(t,N):
for ii in prange(N[0]):
rec_fill(t[ii],N)
#La fonction Ref est un algorithme de type Moulinec et Suquet qui permet de calculer la conductivit homognise (kref) partir d'une 'image' (ki, qui est un tableau contenant les conductivits de chaque pixel de l'image). N est la taille de l'image. tau_field est le champ de polarisation qui sert de point de dpart de l'algorithme. Et k0 est la conductivit de rfrence utilise pour le calcul.
#La fonction Ref est un algorithme de type Moulinec et Suquet qui permet de calculer la conductivit homognise (kref) partir d'une 'image' (ki, qui est un tableau contenant les conductivits de chaque pixel de l'image). N est la taille de l'image. tau_field est le champ de polarisation qui sert de point de dpart de l'algorithme. Et k0 est la conductivit de rfrence utilise pour le calcul.
def Ref(tau_field,N,k0,ki):
d=len(N)
filter_level=2
Eps_f=np.zeros(tuple(N)+(d,d))
prevEpsf=np.zeros(tuple(N)+(d,d))
for j in range(d):
tau = np.zeros(tuple(N)+(d,))
for i in np.ndindex(tuple(N)):
......@@ -50,7 +51,7 @@ def Ref(tau_field,N,k0,ki):
for i in np.ndindex(tuple(N)):
Eps_field[i]+=field[i]
for i in np.ndindex(tuple(N)):
tau[i]=np.dot(ki[i]-k0,Eps_field[i])
tau[i]=np.dot(ki[i]-k0, Eps_field[i])
it+=1
preverr=err
prevkrefjj=krefj[j]
......@@ -61,7 +62,6 @@ def Ref(tau_field,N,k0,ki):
print(krefj)
for i in np.ndindex(tuple(N)+(d,)):
Eps_f[i+(j,)]=Eps_field[i]
prevEpsf[i+(j,)]=prevEps[i]
for i in np.ndindex(tuple(N)):
Sig_field[i]=np.dot(ki[i],Eps_field[i])
plt.figure()
......
{
"cells": [
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"ERROR:root:\n",
"UnicodeDecodeError while processing traceback.\n",
"\n"
]
},
{
"ename": "ImportError",
"evalue": "cannot import name 'jitclass' from 'numba' (/home/bleyerj/.local/lib/python3.8/site-packages/numba/__init__.py)",
"output_type": "error",
"traceback": [
"\u001b[0;31m---------------------------------------------------------------------------\u001b[0m"
]
}
],
"source": [
"from microstructure import booleanEllipses\n",
"from Moulinec_Suquet import Ref\n",
"import numpy as np\n",
"\n",
"d = 2 # dimension (2 components of gradient in 2D)\n",
"N0 = 10 \n",
"N = np.array([N0, N0]) # size of the image (N0 x N0)\n",
"\n",
"centers = [[N0/2, N0/2]] # list of center coordinates in image units\n",
"radius = [N0/4] # list of radius for circles\n",
"theta = [0] # list of orientations (used only for ellipses)\n",
"microstructure = booleanEllipses(centers, radius, radius, theta, N)\n",
"print(microstructure.phasemap)\n",
"\n",
"\n",
"k1 = 1 # conductivity of phase 0\n",
"k2 = 10 # conductivity of phase 1\n",
"k0 = (k1+k2)/2*np.eye(d) # conductivity of reference medium\n",
"\n",
"# initialize conductivity field\n",
"ki=np.zeros(tuple(N) + (d, d))\n",
"# initialize polarization field\n",
"tau=np.zeros(tuple(N) + (d, d))\n",
"for i in np.ndindex(tuple(N)):\n",
" chi = microstructure.phasemap[i] # characteristic function of pixel i\n",
" ki[i] = ((1-chi)*k1 + chi*k2)*np.eye(d)\n",
" tau[i]=np.dot(ki[i]-k0,np.eye(d))\n",
" \n",
"\n",
"kref=Ref(tau,N,k0,ki)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.8.5"
}
},
"nbformat": 4,
"nbformat_minor": 5
}
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment