Projets

Modélisation de l’îlot de chaleur urbain

Ce début d’année 2021 voit l’aboutissement du développement d’un outil de simulation des températures de surfaces en milieu urbain par une approche en couplage aéraulique faible.

L’outil permet le calcul des indicateurs de confort adaptés à l’environnement extérieur en chaque point de l’espace. Il intègre le couplage de simulations de mécanique des fluides, d’ensoleillement et de modèles de parois en régime dynamique.

Les résultats d’un cas d’étude sur le quartier de la gare de Strasbourg sont présentés ci-après.

L’animation suivante montre l’évolution des niveaux de confort sur la première semaine de l’année. Sit back and relax!

Animation GIF Strasbourg ICU
Évolution du niveau de confort sur une semaine de janvier.
blog, Publication

Publication dans les Techniques de l’Ingénieur

Publication de notre méthode de calcul pour la ventilation naturelle

Notre approche de couplage entre mécanique des fluides numériques et simulation thermique dynamique, servant au calcul de ventilation naturelle, a été publiée ce mois-ci dans les Techniques de l’Ingénieur. L’article complet est à découvrir sur leur site. Alternativement vous pouvez nous contacter pour plus de précisions !

Projets

Gare de Paris Est

Variantes de couverture et d’aménagement : influence sur le confort thermique intérieur

L’étude du niveau de confort dans le projet d’aménagement de la gare de l’Est a été réalisée en utilisant plusieurs méthodes :

  • Un couplage entre CFD et simulation thermique dynamique pour le calcul de la ventilation naturelle avec des coefficients de pression calculés finement
  • La distribution intérieure détaillée des flux solaires avec le logiciel Radiance
  • L’obtention des niveaux de confort en intérieur heure par heure pour toute l’année puis leur traitement statistique

Quelques résultats graphiques ci-dessous.

/

Projets

Gare de Lyon Part Dieu

Étude de la distribution intérieure du confort thermique

L’étude du niveau de confort a été réalisée en utilisant plusieurs méthodes :

  • Un couplage entre CFD et simulation thermique dynamique pour le calcul de la ventilation naturelle avec des coefficients de pression calculés finement
  • La distribution intérieure détaillée des flux solaires avec le logiciel Radiance
  • L’obtention des niveaux de confort en intérieur heure par heure pour toute l’année à l’aide d’un métamodèle afin d’accélérer les calculs sur ce modèle très volumineux

Le résultat en images !

Perspective du projet (vue d’architecte).

Aperçu du site dans son environnement urbain, utilisé pour le calcul CFD.

Météorologie du site de Lyon-Part-Dieu prenant en compte l’effet d’îlot de chaleur urbain (ICU)

Pourcentage du temps en zone de confort (indicateur OTCA).

Niveau de confort P.E.T. en moyenne annuelle.

Distribution du champ des vitesses moyen en intérieur.

Distribution moyenne de la température radiante.

Projets

Gare de Paris Austerlitz

Modélisation annuelle de confort spatial en intérieur en Gare d’Austerlitz

Une étude sur les variantes possibles pour l’amélioration du  confort d’été et d’hiver du halls principal mettant en œuvre  les phases suivantes :

  • Un couplage CFD/simulation thermique dynamique pour le calcul de la ventilation naturelle avec les bons coefficients de pression
  • Le calcul détaillé des flux solaires en intérieur avec Radiance
  • La détermination des niveaux de confort en intérieur heure par heure pour toute l’année

Projets

Gare de l’Est

Étude annuelle de confort spatial en intérieur

Une étude sur les variantes possibles pour l’amélioration du  confort d’été et d’hiver des halls latéraux ainsi que du quai transversal mettant en œuvre :

  1. Un couplage CFD/STD et flux solaires pour le calcul spatial des niveaux de confort en intérieur.
  2. Une analyse de sensiblité selon la méthode de Morris sur objectif de confort pour déterminer les leviers d’action les plus importants
  3. Le recours à l’optimisation génétique afin d’obtenir les jeux de paramètres maximisant les conforts d’été et d’hiver.

Ci-dessous une animation sur l’évolution des objectifs de confort lors de l’optimisation génétique, en fonction de l’avancée des générations :

Optimisation génétique : avancée des solutions optimales au fil des générations.

 

Ainsi, chaque « point » qui apparait correspond à une étude complète de confort (point 1 ci-dessus) sur le projet. L’objectif de l’optimisation génétique est de trouver le meilleur compromis entre le confort d’été et le confort d’hiver (front de Pareto).

Non classé

Station de mesure autonome

Une station de mesure autonome pour le diagnostic

Fruit d’un partenariat de plus de deux ans entre AREP la spécialité Génie Electrique de l’INSA de Strasbourg, voici une station de mesure autonome, démarrée lors de projets antérieurs (voir dans les références de l’Hypercube) qui a abouti en janvier 2019 grâce au Projet de Recherche Technologique de François-Alexandre Fournier, élève ingénieur en GE5 à l’INSA : nous le remercions ici pour son investissement sans égal !

Comportant une base Raspberry PI et des capteurs basse consommation sans fil Whisper node de WISEN, cette station développée sur mesure pour les besoins d’AREP comporte notamment une dizaine de sondes de températures et d’humidité ainsi qu’un capteur de CO2. Un travail particulier sur la consommation d’énergie et la robustesse de fonctionnement a été mené au cours du projet. Essais de terrain prévus en 2019.

[Station autonome : Raspberry PI + communication 868 MHz avec les capteurs sans fil]

Projets

Gare de St-Michel Notre Dame (RER C)

Développement d’un modèle physique à équations différentielles pour la prédiction de la concentration en PM10 dans les gares souterraines, dans le but d’évaluer la pertinence de solutions techniques d’amélioration de la qualité d’air (ventilation, filtration, …).

Pour plus de détails voir les publications issues des travaux sur ce modèle (article dans TRD, article dans STOTEN).

Non classé

Projection de champ

Une routine en python pour projeter deux champs sur une grille régulière en 2D, selon une adaptation simplifiée de l’algorithme des plus proches voisins d’Aster, aussi détaillé sur notre page. Le code est paramétré en nombre de points alentours et exécutable en parallèle, par exemple en utilisant multiprocessing de python.

##################################################################
# _______  ______    _______  _______ 
#|   _   ||    _ |  |       ||       |
#|  |_|  ||   | ||  |    ___||    _  |
#|       ||   |_||_ |   |___ |   |_| |
#|       ||    __  ||    ___||    ___|
#|   _   ||   |  | ||   |___ |   |    
#|__| |__||___|  |_||_______||___|    
#
##################################################################
# ce sctript permet la projection de deux champs sur une grille  #
# régulière                                                      #
##################################################################
# Last modification :  14/08/2018                                #
##################################################################
# Copyright (C) 2018 Edouard WALTHER                          #  
# This program is free software; you can redistribute it and/or  #
# modify it under the terms of the GNU General Public License    #
# as published by the Free Software Foundation; either version   #
# of the License, or (at your option) any later version.         #
##################################################################
# Contact : edouard[dot]walther[at]arep[dot]com                                
##################################################################


import numpy as np
import math

def calc_parallel(argu_heure):

	print "quelle heure", argu_heure
	# fichier source
	dossier_src="./maillage/"
	nom_fichier=dossier_src+"coordonnees.txt"

	# ouverture
	coord=open(nom_fichier,"r")
	lignes=coord.readlines()
	coord.close()
	# on les met dans trois vecteurs uniques
	x=[]
	y=[]

	# comme les x,y sont en colonnes on les separe
	for line in lignes:
		columns = line.split(",") # separateur = virgule+espace(s)
		if not columns[0].startswith("#"): # verif inutile mais c'est cadeau
			x.append(float(columns[0]))
			y.append(float(columns[1]))

	# le type de flukx
	typeflux="diffus"
	dossier_src="./flux_"+typeflux+"/"
	dossier_dst="./flux_"+typeflux+"_projete/"
	suffixe="simulation_annuelle_"+typeflux+"_"
			
	# on definit les min / max des coordonnees
	xmin = min(x)
	xmax = max(x)
	ymin = min(y)
	ymax = max(y)
	#zmin = min(z)
	#zmax = max(z)

	#nb reel d'elements
	N_elts=len(x)
	#dx_reel_moyen=(xmax-xmin)/float(Nx_reel)
	
	# nb d'elements pour la grille reguliere
	Nx=180
	Ny=180

	# les pas reguliers d'espace (dx,dy)
	dx=(xmax-xmin)/float(Nx-1)
	dy=(ymax-ymin)/float(Ny-1)

	# print pour pouvoir en faire qque chose
	# print "Les coordonnees"
	# print "	xmin, xmax", xmin, xmax
	# print "	ymin, ymax", ymin, ymax
	# print "	dx,dy", dx, dy

	# les matrices contenant x,y reguliers (on peut faire sans en recodant)
	global grille_x
	global grille_y
	grille_x=np.ones((Nx,Ny))
	grille_y=np.ones((Nx,Ny))
	# une matrice pour verifier qu'on a bien rempli la grille
	global grille_verif
	grille_verif=np.zeros((Nx,Ny))

	# je mets ici l'exposant pour le lissage
	beta=0.333 # 2*0.75
	# le nb de points voisins pour voir
	nb_pts=4
	# une fc pour calculer la distance entre deux points 
	# (ca soulage les yeux dans la partie d'apres!)
	# il n'y a que 3 coordonnees car x et y ont le meme indice "i"
	def dist(i,x2,y2):
		dist_xy=math.sqrt(math.pow(x[i]-grille_x[x2][y2],2)+math.pow(y[i]-grille_y[x2][y2],2) )
		return dist_xy

	# creation grille reguliere
	for k in range(Nx):
		for l in range(Ny):
			grille_x[k][l] = k*dx+xmin
			grille_y[k][l] = l*dy+ymin
			
	for heure in range(argu_heure,argu_heure+1):
		# ponderations
		ponderation=np.zeros((Nx,Ny))
		somme_dist=np.zeros((Nx,Ny))
		champ=np.ones((Nx,Ny))*999.99
		
		print " le fichier", suffixe+str(heure)+".txt"
		
		nom_fichier_source=dossier_src+suffixe+str(heure)+".txt"
		nom_fichier_destination=dossier_dst+suffixe+str(heure)+".txt"
		fichier = open(nom_fichier_source, "r")
		z=map(float,fichier) #pour float > string
		fichier.close()
		# si z ne contient rien (pas de flux, genre la nuit)
		if z==[]:
			print "	(est vide)"
			fichier=open(nom_fichier_destination,"w")
			fichier.close()
		else:
			print " (n'est pas vide)"
			# balayer les points du nuage
			for i in range(N_elts):		
				# on trouve le point du maillage regulier le plus proche 
				# 	>> sur l'axe x
				if (x[i]- xmin) % dx > dx/2.:
					i_grille=1+int(math.floor( (x[i]-xmin)/dx ))
				else:
					i_grille=int(math.floor( (x[i]-xmin)/dx ))
				# 	>> sur l'axe y
				if (y[i]- ymin) % dy > dy/2.:
					j_grille=1+int(math.floor( (y[i]-ymin)/dy ))
				else:
					j_grille=int(math.floor( (y[i]-ymin)/dy ))

				# on stocke les valeurs de la grille qui ont ete remplies
				grille_verif[i_grille][j_grille]=1. # verif que touchee

				for m in range(-nb_pts,nb_pts+1):
			        for n in range(-nb_pts,nb_pts+1):
				        if i_grille+m>0 and j_grille+n>0 and j_grille+n<Ny and i_grille+m<Nx and grille_touched[i_grille+m][j_grille+n]==1:
        					distance=dist(i,i_grille+m,j_grille+n)
        					weight=math.exp(-math.pow(distance/dref,beta))
        					if beta==0:weight=1
        					ponderation[i_grille+m][j_grille+n]=ponderation[i_grille+m][j_grille+n]+weight*grille_touched[i_grille+m][j_grille+n]*z[i]
        					somme_dist[i_grille+m][j_grille+n]=  somme_dist[i_grille+m][j_grille+n]+weight*grille_touched[i_grille+m][j_grille+n]
        					grille_verif[i_grille+m][j_grille+n]=1

			print "Stockage dans la matrice"
			for i in range(Nx):
				for j in range(Ny):
					if grille_verif[i][j]==0: # si aucun point n'a ete touche
						zob=1
						#champ[i][j]=-5555.
					else:
						champ[i][j]= ponderation[i][j]/somme_dist[i][j]
			fichier=open(dossier_dst+suffixe+ str(heure)+".txt","w") 
			for i in range(Nx):
				for j in range(Ny):
					chaine=str(grille_x[i][j])+" "+str(grille_y[i][j])+" "+str(champ[i][j])
					fichier.write(chaine+"\n") # retour ligne
				fichier.write("\n") # saut de ligne entre les x
			fichier.close()<span id="mce_marker" data-mce-type="bookmark" data-mce-fragment="1">​</span>

Une deuxième version plus aboutie consiste à résoudre un problème de minimisation pour projeter le champ sur une grille (ici régulière) à partir des valeurs du nuage de points. Le code proposé ci-dessous se sert de la librairie scipy.optimize.

# _______  ______    _______  _______ 
#|   _   ||    _ |  |       ||       |
#|  |_|  ||   | ||  |    ___||    _  |
#|       ||   |_||_ |   |___ |   |_| |
#|       ||    __  ||    ___||    ___|
#|   _   ||   |  | ||   |___ |   |    
#|__| |__||___|  |_||_______||___|    
#
##################################################################
# ce sctript permet la projection de deux champs sur une grille  #
# régulière                                                      #
##################################################################
# Last modification :  16/01/2019                                #
##################################################################
# Copyright (C) 2018 Edouard WALTHER   / Antoine HUBERT          #  
# This program is free software; you can redistribute it and/or  #
# modify it under the terms of the GNU General Public License    #
# as published by the Free Software Foundation; either version   #
# of the License, or (at your option) any later version.         #
##################################################################
# Contact : edouard[dot]walther[at]arep[dot]com                                
##################################################################

import os
import numpy as np
import time
import math
import csv
from scipy.optimize import minimize
debut=time.time()
alsace=True

algo="CG"
algo="SLSQP"
direction=1
algo="POWELL"
# nb de points en x et y
Nx=30
Ny=30
# la grille touchee (une matrice Nx*Ny remplie de 1 fait l'affaire si pas de "trous" dans votre maillage)
grille_touched=np.load("grille_touched_"+str(Nx)+ '.npy')
les_min_max=np.load("min_max"+ '.npy') #les min max des coordonnees
minxparoi1=les_min_max[0]
maxxparoi1=les_min_max[1]
minyparoi1=les_min_max[2]
maxyparoi1=les_min_max[3]

# quelques defs initiales pour les dossiers source et destination
if alsace==True:prefixe_dst="./vitesse_existant_hallHA"
else:prefixe_dst="./vitesse_existant_hallSM"

argu_heure=0
dossier_src="../../RESULTATS/hall_alsace/CSV_point/"
dossier_dst="./minimisation_matrice_beta_dref_flux_"+algo+"/"
#creation
if not os.path.isdir(dossier_dst):
	os.mkdir(dossier_dst)
#lecture des fichiers
if alsace==True:
	fichier_source=dossier_src+'flux_point.'+str(argu_heure)+".csv"

# def des coordonnees du champ : on prend un peu de marge
xmin = minxparoi1-4
xmax = maxxparoi1+4
ymin = minyparoi1-4
ymax = maxyparoi1+4
dx=(xmax-xmin)/float(Nx-1)
dy=(ymax-ymin)/float(Ny-1)
grille_x=np.ones((Nx,Ny))
grille_y=np.ones((Nx,Ny))

#une fc de calcul de la distance qui soulage les yeux
def dist(i,x2,y2):
	return math.sqrt(math.pow(x[i]-grille_x[x2][y2],2)+math.pow(y[i]-grille_y[x2][y2],2) )

# creation grille reguliere
for k in range(Nx):
	for l in range(Ny):
		grille_x[k][l] = k*dx+xmin
		grille_y[k][l] = l*dy+ymin
# lecture du fichiers
x,y,z = [],[],[]
with open(fichier_source,'r') as csvfile:
		reader = csv.reader(csvfile,delimiter=',',quoting = csv.QUOTE_NONNUMERIC)
		next(reader)
		for line in reader:
			x.append(line[1])
			y.append(line[2])
			z.append(line[0])

N_elts=len(x)
#-----------------------------------------------------
#
#	la fc de projection selon les plus proches voisins comme defini dans code Aster
#
def fc_projection_deg_1(vecteur_params):
	
	a=vecteur_params[0:Nx*Ny]
	b=vecteur_params[Nx*Ny:2*Nx*Ny]
	c=vecteur_params[2*Nx*Ny:3*Nx*Ny]
	beta=1.5 # exposant de "lissage" de la ponderation
	numerateur_dref=3.# le numerateur de dref=dx/numerateur_dref
	nb_pts=4 # cbien de points autour va-t-on chercher
	# ponderations
	dref=dx/numerateur_dref
	# balayer les points du nuage -----------------
	somme_moindre_ca=0.
	for i in range(N_elts):
	    # le carre le plus proche
		i_grille=int(math.floor((x[i]-xmin)/dx))
		j_grille=int(math.floor((y[i]-ymin)/dy))
		# on regarde les points autour
		for m in range(-nb_pts,nb_pts+1):
			for n in range(-nb_pts,nb_pts+1):
				if i_grille+m>0 and j_grille+n>0 and j_grille+n<Ny and i_grille+m<Nx and grille_touched[i_grille+m][j_grille+n]==1:
					distance=dist(i,i_grille+m,j_grille+n)
					weight=math.exp(-math.pow(distance/dref,beta))
					b_x = b[(i_grille+m)*Ny+j_grille+n]*grille_x[i_grille+m][j_grille+n]
					c_y = c[(i_grille+m)*Ny+j_grille+n]*grille_y[i_grille+m][j_grille+n]
					F= a[(i_grille+m)*Ny+j_grille+n] + b_x + c_y
					somme_moindre_ca=somme_moindre_ca+ weight*(F-z[i])**2
	return somme_moindre_ca

# les defs du 
a_tot=np.ones((3*Nx*Ny+2))
a_tot[-2]=1.5 #beta
a_tot[-1]=3. #numerateur_dref

# si on a une estimation precedente on peut la charger
#vec_precedent="./minimisation_matrice_beta_dref_flux_POWELL/vecteur_abc.npy"
#a_tot=np.load(vec_precedent)

# def des limites pour pour les coeffs a_k,b_k,c_k
#	et pour beta, num_dref
bnds=()
for k in range(3*Nx*Ny):
	bnds=bnds+((None,None),)
#limites pour beta 
bnds=bnds+((1,5),)
# ...et pour numerateur
bnds=bnds+((1,10),)

print("Debut minimisation pour l'algo " +algo + "...")

if algo=="POWELL":
	res=minimize(fc_projection_deg_1, a_tot, method='Powell', tol=1e-5, options={'disp':True, 'xtol':0.0001,'ftol':0.0001,'maxiter':10,'maxfev':30000})
if algo=="SLSQP":
	res=minimize(fc_projection_deg_1, a_tot, method='SLSQP', tol=1e-6, bounds=bnds, options={'disp':True,'ftol':0.00001,'maxiter':90})
if algo=="CG":
	res=minimize(fc_projection_deg_1, a_tot, method='CG', tol=1e-6, bounds=bnds, options={'disp':True, 'gtol':0.0001,'maxiter':10})
print(res.message)
print("Vecteur x", res.x)

#on sauvegarde le resultat
matrice=dossier_dst+"vecteur_abc.npy"
np.save(matrice, res.x)

#ensuite on ecrit sous le bon format pour gnuplot
vecteur_params=np.load(matrice)
a=vecteur_params[0:Nx*Ny]
b=vecteur_params[Nx*Ny:2*Nx*Ny]
c=vecteur_params[2*Nx*Ny:3*Nx*Ny]
print ("Stockage...")
fichier=open(dossier_dst+'flux_'+str(argu_heure)+".txt","w") 
for i in range(Nx):
	for j in range(Ny):
		# point touche ET dans la zone d'interet originelle
		if grille_touched[i][j]==0:
			chaine=str(grille_x[i][j])+" "+str(grille_y[i][j])+" "+str(999.99)
		else:
			b_x = b[i*Ny+j]*grille_x[i][j]
			c_y = c[i*Ny+j]*grille_y[i][j]
			F= a[i*Ny+j] + b_x + c_y
			chaine=str(grille_x[i][j])+" "+str(grille_y[i][j])+" "+str(F)
		fichier.write(str(chaine)+"\n") # retour ligne
	fichier.write("\n") # saut de ligne entre les x
fichier.close()
fin=time.time()

# calcul des temps passes
deltaT=round(fin-debut,1)
suffixe=" (s)"
if deltaT > 60:
	deltaT=round(deltaT/60,1)
	suffixe=" (min)"
if deltaT > 60:
	deltaT=round(deltaT/60,1)
	suffixe=" (h)"
print ("le calcul a duré "+str(deltaT)+suffixe +" pour l'algorithme " + algo)