On the numerical solution of compressible fluid flow and radiative heat transfer problems

Abstract

Σκοπός της παρούσας Διδακτορικής Διατριβής ήταν η ανάπτυξη μεθοδολογίας για την αριθμητική επίλυση προβλημάτων μόνιμης ροής συμπιεστού ρευστού και μετάδοσης θερμότητας μέσω ακτινοβολίας. Η υπόψη μέθοδος δύναται να χρησιμοποιηθεί σε ένα ευρύ φάσμα εφαρμογών τόσο σε ακαδημαϊκό όσο και σε βιομηχανικό επίπεδο, π.χ., σε θαλάμους καύσης, βιομηχανικούς φούρνους, εσωτερική ή εξωτερική αεροδυναμική ροή, κινητήρες τύπου ramjet, κ.λπ. Πιο συγκεκριμένα, αναπτύχθηκε μεθοδολογία τύπου χρονοπροέλασης για την προσομοίωση ατριβούς και συνεκτικής στρωτής ή τυρβώδους μόνιμης συμπιεστής ροής. Η διακριτοποίηση των εξισώσεων Navier-Stokes επί τρισδιάστατων τετραεδρικών ή υβριδικών μη-δομημένων υπολογιστικών πλεγμάτων επιτυγχάνεται με την κεντροκομβική μέθοδο πεπερασμένων όγκων (FVM/Finite Volume Method). Για την μοντελοποίηση τυρβώδους ροής εφαρμόζονται οι κατά Reynolds ολοκληρωμένες εξισώσεις Navier-Stokes (RANS/Reynolds-Averaged Navier-Stokes), χρησιμοποιώντας την υπόθεση Boussinesq και επακόλουθα τον όρο της τυρβώδους συνεκτικότητας, για τον υπολογισμό της οποίας συμπεριλήφθηκαν τρία διαφορετικά μοντέλα τύρβης, το k-ε, το k-ω και το SST (Shear Stress Transport). Για την εκτίμηση των μη-συνεκτικών διανυσμάτων ροής εφαρμόζεται ο προσεγγιστικός επιλύτης του Roe, θεωρώντας ένα τοπικό μονοδιάστατο πρόβλημα Riemann στη διεπαφή των γειτονικών όγκων ελέγχου. Αύξηση της ακρίβειας του προαναφερθέντος υπολογισμού επιτυγχάνεται με την εφαρμογή σχήματος δεύτερης τάξης χωρικής ακρίβειας, βασισμένου στην τεχνική MUSCL (Monotonic Upwind Scheme for Conservation Laws). Το εν λόγω σχήμα συνδυάζεται με κατάλληλη συνάρτηση περιορισμού (Van Albada-Van Leer, Min-mod ή Barth-Jespersen) προκειμένου να διασφαλιστεί η μονοτονία μεταξύ των τιμών των μεταβλητών των γειτονικών όγκων ελέγχου. Η εκτίμηση των συνεκτικών διανυσμάτων ροής προϋποθέτει τον πρωτύτερο υπολογισμό των παραγώγων των συνιστωσών της ταχύτητας και της θερμοκρασίας στη διεπαφή των όγκων ελέγχου, η οποία συμπίπτει με το μέσο της ακμής που συνδέει τους αντίστοιχους υπολογιστικούς κόμβους. Για τον υπόψη υπολογισμό εισήχθησαν δύο τεχνικές στην παρούσα μεθοδολογία, εκ των οποίων η πρώτη βασίζεται στη δημιουργία νέων δυικών όγκων ελέγχου γύρω από την υπό εξέταση ακμή (edge-dual volume method), ενώ σύμφωνα με τη δεύτερη μέθοδο οι επιθυμητές παράγωγοι προκύπτουν από τις αντίστοιχες των ακραίων κόμβων της ακμής (nodal-averaging method). Ο υπολογισμός τόσο των μη-συνεκτικών όσο και των συνεκτικών διανυσμάτων ροής εκτελείται με σάρωση των ακμών του πλέγματος, χρησιμοποιώντας κατάλληλες δομές δεδομένων (edge-based data structures), προκειμένου να μειωθεί όσο το δυνατόν ο απαιτούμενος υπολογιστικός χρόνος. Η χρονική ολοκλήρωση και τελική κατάσταση της ροής προσεγγίζεται επαναληπτικά, είτε με ρητό σχήμα, εφαρμόζοντας την μέθοδο Runge-Kutta τεσσάρων βημάτων (RK(4)) και δεύτερης τάξης χρονικής ακρίβειας, είτε με σημειακά πεπλέγμενο σχήμα, εφαρμόζοντας τον αλγόριθμο Jacobi ή τον αλγόριθμο Gauss-Seidel. Για την επιτάχυνση της επίλυσης εφαρμόζεται επιπρόσθετα η τεχνική του τοπικού ψευδο-χρονικού βήματος (local time-stepping technique). Τέλος, σημειώνεται ότι για τις εξισώσεις των μοντέλων τύρβης ακολουθείται παρόμοια με τις εξισώσεις ροής διαδικασία χρονικής ολοκλήρωσης και υπολογισμού των διανυσμάτων ροής, εκτός του προσεγγιστικού επιλύτη του Roe και του σχήματος δεύτερης τάξης χωρικής ακρίβειας, καθώς η κύρια συνεισφορά στα εν λόγω μοντέλα προέρχεται από τους συνεκτικούς όρους.Για την μοντελοποίηση της μετάδοσης θερμότητας μέσω ακτινοβολίας σε γκρι μέσο με δυνατότητα απορρόφησης, εκπομπής και σκέδασης εφαρμόζεται η εξίσωση μεταφοράς ακτινοβολίας (RTE/Radiative Transfer Equation) με παρουσία (time-dependent) ή μη (non time-dependent) χρονικού όρου. Καθώς η υπόψη εξίσωση θα πρέπει να επιλυθεί για κάθε πεπερασμένο όγκο ελέγχου και για κάθε πεπερασμένη στερεά γωνία ελέγχου, διακριτοποιείται χωρικά και χρονικά. Η εν λόγω χωρική διακριτοποίηση επιτυγχάνεται με τη κεντροκομβική μέθοδο πεπερασμένων όγκων, ομοίως των εξισώσεων ροής, ενώ η γωνιακή παρομοιάζεται με το χωρισμό της επιφάνειας «σφαίρας» (περιμέτρου «κύκλου» στις δύο διαστάσεις) σε μικρότερες αντίστοιχες των επιθυμητών στερεών γωνιών ελέγχου. Για τον υπολογισμό των διανυσμάτων θερμορροής στη διεπαφή των όγκων ελέγχου εφαρμόζεται απλό άναντες σχήμα (step scheme), σύμφωνα με το οποίο η ένταση της ακτινοβολίας στη διεπαφή τίθεται ίση με αυτή του κόμβου στα ανάντη της θερμορροής. Για την αύξηση της ακρίβειας της λύσης χρησιμοποιείται σχήμα δεύτερης τάξης χωρικής ακρίβειας, αντίστοιχο αυτού για τα μη-συνεκτικά διανύσματα ροής των εξισώσεων Navier-Stokes. Το υπόψη σχήμα συνδυάζεται με τη συνάρτηση περιορισμού Van Albada-Van Leer ή Min-mod, προκειμένου να διασφαλιστεί η μονοτονία μεταξύ των τιμών της έντασης της ακτινοβολίας των γειτονικών κόμβων, ιδιαίτερα κοντά στις οριακές επιφάνειες. Το πρόβλημα επικάλυψης όγκων ελέγχου και στερεών γωνιών ελέγχου, στο οποίο αναπόφευκτα οδηγείται η εν λόγω μέθοδος λόγω συνδυασμού μη-δομημένων πλεγμάτων και γωνιακής διακριτοποίησης, αντιμετωπίζεται με την εφαρμογή της Bold Approximation ή της Pixelation Technique. Αν και κατά τη συνήθη τακτική οι οριακές συνθήκες επιβάλλονται ρητά σε τέτοιους είδους μεθοδολογίες, στην παρούσα εργασία υιοθετήθηκε η πεπλεγμένη επιβολή τους, επιτρέποντας τη χρήση αραιότερων πλεγμάτων στην περιοχή των οριακών επιφανειών, καθώς και την ανάπτυξη και εφαρμογή οριακών συνθηκών τύπου συμμετρίας. Τέλος, η προσέγγιση της τελική λύσης της εξίσωσης επιτυγχάνεται είτε με επαναληπτική διόρθωση των τιμών (non time-dependent RTE) είτε με χρονική ολοκλήρωση και εφαρμογή της μεθόδου Runge-Kutta τεσσάρων βημάτων (RK(4)) και δεύτερης τάξης χρονικής ακρίβειας (time-dependent RTE). Με στόχο την επιτάχυνση της αριθμητικής επίλυσης τόσο των προβλημάτων ροής όσο και των προβλημάτων μετάδοσης θερμότητας μέσω ακτινοβολίας (επιπλέον της διάρθρωσης των δεδομένων κατά τις ακμές του υπολογιστικού πλέγματος και της εφαρμογής της τεχνικής του τοπικού ψευδο-χρονικού βήματος), χρησιμοποιήθηκε μεθοδολογία παράλληλης επεξεργασίας και πολυπλέγματος. Η μέθοδος παραλληλοποίησης βασίζεται στην τεχνική διαμέρισης του υπολογιστικού πλέγματος σε μικρότερα υποχωρία (domain decomposition approach), ώστε να καταστεί δυνατή η ταυτόχρονη επίλυση των εξισώσεων σε αυτά. Η όλη διαδικασία ξεκινάει στον κύριο επεξεργαστή με την εφαρμογή του λογισμικού METIS, το οποίο διαχωρίζει τους κόμβους του αρχικού πλέγματος σε μικρότερα σύνολα κόμβων (core nodes). Ωστόσο, με αυτόν τον τρόπο τα στοιχεία του πλέγματος στα όρια των υπο-συνόλων (υπο-πλεγμάτων) παραμένουν ανολοκλήρωτα, καθώς δεν περιέχονται όλοι οι κόμβοι τους στα αντίστοιχα σύνολα. Για την ολοκλήρωση τους προστίθενται οι ελλείποντες κόμβοι (ghost nodes) στα αντίστοιχα υπο-πλέγματα, δημιουργώντας ταυτόχρονα μία περιοχή αλληλοκάλυψης ανάμεσα τους, ενώ στη συνέχεια διανέμονται τα απαραίτητα δεδομένα και στους υπόλοιπους επεξεργαστές. Η ανταλλαγή πληροφορίας μεταξύ των γειτονικών υπο-πλεγμάτων, που αφορά στις τιμές των μεταβλητών καθώς και των παραγώγων αυτών, επιτυγχάνεται μέσω των κόμβων στις περιοχές αλληλοκάλυψης και τις εντολές του πρωτόκολλου επικοινωνίας MPI (Message Passing Interface). Οι προαναφερθείσες τιμές των ghost κόμβων δεν υπολογίζονται από τις εξισώσεις μεταφοράς, αλλά λαμβάνονται απευθείας από τους αντίστοιχους core κόμβους των γειτονικών υπο-πλεγμάτων.Η μέθοδος πολυπλέγματος (multigrid method) στην παρούσα Διδακτορική Διατριβή εφαρμόζεται σε χωρική μορφή για την περίπτωση προβλημάτων ροής και μετάδοσης θερμότητας μέσω ακτινοβολίας, καθώς και σε γωνιακή και συνδυασμένη χωρική-γωνιακή μορφή μόνο για προβλήματα ακτινοβολίας. Η χωρική μέθοδος βασίζεται στη δημιουργία μίας σειράς αραιότερων του αρχικού πλεγμάτων και στην επίλυση προσεγγιστικών εξισώσεων σε αυτά, με σκοπό την ταχύτερη σύγκλιση στην τελική λύση. Για τη δημιουργία αυτών των πλεγμάτων εφαρμόζεται η τεχνική συσσωμάτωσης (agglomeration approach), σύμφωνα με την οποία οι όγκοι ελέγχου των γειτονικών κόμβων ενώνονται δημιουργώντας ένα νέο υπερ-κόμβο του αραιότερου υπο-πλέγματος. Η υπόψη συσσωμάτωση, η οποία δύναται να είναι είτε ισότροπη (isotropic) είτε κατευθυνόμενη (directional), πραγματοποιείται με τρόπο τέτοιο ώστε να διατηρεί την αρχική τοπολογία, ενώ ξεκινάει από τις οριακές επιφάνειες και εκτείνεται προς το εσωτερικό του πλέγματος προσομοιάζοντας την τεχνική του προελαύνοντος μετώπου (advancing front technique). Ωστόσο, η εν λόγω ένωση υπόκειται σε προκαθορισμένους περιορισμούς, που αφορούν κυρίως στους κόμβους των εσωτερικών και εξωτερικών ορίων των υπο-πλεγμάτων, π.χ., ένας οριακός κόμβος δύναται να ενωθεί μόνο με έναν άλλο οριακό του ίδιου είδους, ενώ οι "ghost" κομβοι δε λαμβάνονται υπόψη κατά την κύρια διαδικασία αλλά ενώνονται σύμφωνα με τη συσσωμάτωση που έχουν υποστεί οι αντίστοιχοι τους core κόμβοι στα γειτονικά υπο-πλέγματα. Λαμβάνοντας υπόψη τα ανωτέρω, η όλη διαδικασία ξεκινάει σε κάθε υπο-πλέγμα με την ένωση των όγκων ελέγχου των οριακών κόμβων στερεάς επιφάνειας με αυτούς των γειτονικών τους, επίσης οριακών κόμβων, ενώ στη συνέχεια καταρτίζεται λίστα με τους κόμβους που έχουν έρθει σε επαφή με το "μέτωπο" της συσσωμάτωσης (seed nodes). Η ένωση των κόμβων και η δημιουργία υπερ-κόμβων συνεχίζεται με τη συσσωμάτωση των όγκων ελέγχου των "seed" κόμβων με αυτούς των γειτονικών τους. Κατόπιν καταρτίζεται μία νέα λίστα "seed" κόμβων με τον ίδιο τρόπο και η όλη διαδικασία επαναλαμβάνεται έως ότου εξεταστούν όλοι οι κόμβοι του υπολογιστικού πλέγματος. Στην περίπτωση της κατευθυνόμενης συσσωμάτωσης, η οποία εφαρμόζεται σε πλέγματα υβριδικού τύπου και δύναται να είναι είτε μερική (semi-directional) είτε ολική (full-directional), η διαδικασία ξεκινάει από τους οριακούς κόμβους των πρισματικών στοιχείων και συνεχίζεται στους κόμβους των επόμενων πρισματικών στρωμάτων, διατηρώντας το μοτίβο συσσωμάτωσης των πρώτων. Αφού ολοκληρωθεί η ένωση των πρισματικών κόμβων, η διαδικασία συνεχίζεται ισότροπα για τους υπόλοιπους κόμβους, εκκινώντας από τους γειτονικούς των ανώτερων πρισματικών στρωμάτων. Όπως αναφέρθηκε παραπάνω, οι "ghost" κόμβοι δε λαμβάνονται υπόψη κατά την κύρια διαδικασία, αλλά ενώνονται σύμφωνα με τη συσσωμάτωση που έχουν υποστεί οι αντίστοιχοι τους "core" κόμβοι στα γειτονικά υπο-πλέγματα. Με τον τρόπο αυτό δημιουργούνται "ghost" ψευδο-υπερ-κόμβοι, καθώς ο περικλειόμενος αριθμός κόμβων τους πιθανόν να διαφέρει από αυτόν στους αντίστοιχους core υπερ-κόμβους. Εφόσον απαιτείται ακόμη αραιότερο πλέγμα, η όλη διαδικασία επαναλαμβάνεται. Η επίλυση των εξισώσεων στη εν λόγω σειρά πλεγμάτων πραγματοποιείται με το σχήμα FAS (Full Approximation Scheme) σε κύκλο σχήματος V, σύμφωνα με το οποίο στα αραιότερα πλέγματα επιλύεται μία προσεγγιστική μορφή των εξισώσεων μεταφοράς. Η μεταφορά πληροφορίας σε έναν υπερ-κόμβο από τους περικλειόμενους κόμβους του (restriction) αφορά στις σταθμισμένες κατ' όγκο συντηρητικές μεταβλητές τους και στο άθροισμα των διανυσμάτων ροής τους. Αντίθετα, από τον υπερ-κόμβο μεταφέρονται στους περικλειόμενους κόμβους του (prolongation) οι διορθώσεις των μεταβλητών είτε με απλή μεταφορά στην περίπτωση ατριβούς ροής είτε με μεταφορά σταθμισμένη με την απόσταση μεταξύ τους στην περίπτωση συνεκτικής ροής. Για την επίτευξη ακόμη μεγαλύτερης επιτάχυνσης κατά τα πρώτα στάδια του σχήματος FAS, εφαρμόζεται το σχήμα FMG (Full Multigrid).Η γωνιακή μέθοδος πολυπλέγματος εφαρμόζεται με ανάλογο της χωρικής τρόπο στη γωνιακή διακριτοποίηση των προβλημάτων μετάδοσης θερμότητας μέσω ακτινοβολίας. Η συσσωμάτωση συνίσταται στην ένωση των γειτονικών στερεών γωνιών ελέγχου, μειώνοντάς τες στο ένα τέταρτο των αρχικών κάθε φορά. Σε αντίθεση με τη χωρική, η γωνιακή συσσωμάτωση υπόκειται σε έναν μόνο περιορισμό που απαγορεύει την ένωση στερεών γωνιών που ανήκουν σε διαφορετικά τεταρτημόρια. Με αυτόν τον τρόπο διασφαλίζεται η συνέχεια κατά τη μεταφορά της πληροφορίας από την πυκνότερη στην αραιότερη διακριτοποίηση και αντίστροφα. Η επίλυση της RTE στα διαδοχικά αραιότερα γωνιακά "πλέγματα" πραγματοποιείται με το σχήμα FAS ομοίως της χωρικής μεθόδου πολυπλέγματος, ενώ η στάθμιση των μεταβλητών κατά τη μεταφορά τους από το πυκνότερο στο αραιότερο υλοποιείται κατά γωνία (όχι κατ' όγκο). Τέλος, στην παρούσα εργασία αναπτύχθηκε και συνδυασμένη χωρική-γωνιακή μεθοδολογία για την επιτάχυνση προβλημάτων μετάδοσης θερμότητας μέσω ακτινοβολίας. Σύμφωνα με αυτήν, σε κάθε επίπεδο του χωρικού σχήματος πολυπλέγματος εκτελείται ένας πλήρης FAS V-κύκλος του γωνιακού σχήματος, ενώ η μεταφορά της έντασης της ακτινοβολίας, των διανυσμάτων ροής και των διορθώσεων εκτελούνται ομοίως προς τις απλέςν μεθόδους.Επιπρόσθετα, με σκοπό την περαιτέρω αύξηση της ακρίβειας της παρούσας μεθόδου και πέραν της εφαρμογής των προαναφερθέντων σχημάτων υψηλότερης τάξης χωρικής ακρίβειας, αναπτύχτηκε μεθοδολογία αυτόματης τοπικής πύκνωσης του πλέγματος. Με αυτόν τον τρόπο δύνανται να αυξηθούν οι βαθμοί ελευθερίας του υπό εξέταση πλέγματος, αποφεύγοντας ταυτόχρονα την απαίτηση κατασκευής ενός νέου εξαρχής. Η υπόψη τεχνική αναδεικνύεται ιδιαίτερα πολύτιμη σε περιπτώσεις εμφάνισης τοπικών φαινομένων, π.χ. περιοχές κυμάτων κρούσης, αποκόλλησης ροής, κ.λπ., καθώς και σε προβλήματα που αντιμετωπίζονται για πρώτη φορά, με συνέπεια να μην είναι εκ των προτέρων γνωστό το απαιτούμενο επίπεδο πύκνωσης στις διάφορες περιοχές του πλέγματος. Η όλη διαδικασία δύναται να διαιρεθεί σε τέσσερα βασικά βήματα: α) Εντοπισμός των περιοχών προς πύκνωση και σημείωση των αντίστοιχων ακμών βάσει ενός προκαθορισμένου κριτηρίου, π.χ. ακμές που περιλαμβάνουν κόμβους με τιμή του αριθμού Mach υψηλότερη της μονάδας. β) Διάχυση της πληροφορίας σημείωσης στις γειτονικές ακμές, ώστε να καταστεί δυνατή η διαίρεση των αντίστοιχων στοιχείων του πλέγματος. γ) Διαίρεση των σημειωμένων ακμών και εισαγωγή νέων κόμβων στο μέσο τους. δ) Διαίρεση των αντίστοιχων πλευρών και στοιχείων του πλέγματος, βάσει προκαθορισμένων κανόνων διαίρεσης, π.χ., ένα πρισματικό στοιχείο δύναται να διαιρεθεί σε δύο ή τέσσερα νέα πρισματικά στοιχεία.Εφαρμόζοντας τις ανωτέρω τεχνικές αναπτύχθηκε ο ακαδημαϊκός κώδικας Galatea (Γαλάτεια) στα πλαίσια της παρούσας Διδακτορικής Διατριβής. Για την αξιολόγηση του εν λόγω λογισμικού και κατ' επέκταση των προαναφερθεισών μεθόδων εξετάστηκαν διάφορα προβλήματα αναφοράς, ενώ τα αποτελέσματα του συγκρίθηκαν ποιοτικά και ποσοτικά με διαθέσιμα πειραματικά δεδομένα, καθώς και με αριθμητικά αποτελέσματα αντίστοιχων αλγορίθμων αναφοράς. Από τις υπόψη συγκρίσεις διαφαίνεται η δυνατότητα της εν λόγω μεθοδολογίας για τέτοιου είδους προσομοιώσεις τόσο από την πλευρά της ακρίβειας όσο και από την πλευρά της αποδοτικότητας. Εν κατακλείδι, λαμβάνοντας υπόψη τις αντίστοιχες διαθέσιμες στη διεθνή βιβλιογραφία μελέτες, η συνεισφορά της παρούσας εργασίας συνοψίζεται στα κάτωθι:•Στην ανάπτυξη της γωνιακής μεθόδου πολυπλέγματος.•Στην ανάπτυξη της συνδυασμένης χωρικής-γωνιακής μεθόδου πολυπλέγματος.•Στην εφαρμογή σχήματος δεύτερης τάξης χωρικής ακρίβειας με συναρτήσεις περιορισμού για την προσομοίωση της μετάδοσης θερμότητας μέσω ακτινοβολίας.•Στην ανάπτυξη της πεπλεγμένης επιβολής των οριακών συνθηκών σε προβλήματα μετάδοσης θερμότητας μέσω ακτινοβολίας.•Στην τοπική προσαρμογή του πλέγματος για την προσομοίωση της μετάδοσης θερμότητας μέσω ακτινοβολίας.•Στη χρήση υβριδικών πλεγμάτων σε συνδυασμό με την κεντροκομβική μέθοδο πεπερασμένων όγκων σε προβλήματα μετάδοσης θερμότητας μέσω ακτινοβολίας.•Στη χωρική μέθοδο πολυπλέγματος (δημιουργία ghost ψευδο-υπερ-κόμβων, μεταφορά διορθώσεων από το αραιότερο στο πυκνότερο πλέγμα βάσει απόστασης των κόμβων, κ.λπ.).