Σελίδα 1 από 1

Αριθμητικές μέθοδοι

ΔημοσίευσηΔημοσιεύτηκε: 24 Απρ 2010, 13:28
από Dimitris
Μιας και υπάρχουν πολλοί οδηγοί στο forum για προγραμματισμό σε διάφορες γλώσσες, με διάφορες βιβλιοθηκές (Gtk/Qt) και μιας και βαριέμαι για άλλη φορά, είπα να γράψω και έναν οδηγό πιο εξειδικευμένο. Το θέμα αυτού το νήματος είναι οι αριθμητικές μέθοδοι και υλοποίησή τους σε γλώσσες προγραμματισμού. Δεν είναι ένας παραδοσιακός οδηγός (ίσως να πρέπει να αλλαχτεί ενότητα) αλλά έχει ως σκοπό να πληροφορήσει τους ενδιαφερόμενους χρήστες για τις δυνατότητες που υπάρχουν εκεί έξω.

Το πρώτο κεφάλαιο των περισσότερων βιβλίων για αριθμητικές μεθόδους είναι η επίλυση συστημάτων γραμμικών εξισώσεων. Υπάρχουν πάρα πολλές μέθοδοι και δε θα τις σχολιάσω. Το σημαντικό είναι ότι τις περισσότερες φορές η υλοποίηση των μεθόδων, αν δεν γίνεται για εκπαιδευτικούς σκοπούς, υπόκειται σε λάθη και είναι σίγουρα αργή. Για το λόγο αυτό υπάρχει διάφορες βιβλιοθήκες. Μια από αυτές, η πλέον διαδεδομένη, είναι η LAPACK (Linear Algebra PACKage)
http://www.netlib.org/lapack/
Μη βιαστείτε να τη μεταγλωττίσετε. Υπάρχει στα αποθετήρια. Κι ακόμη καλύτερα ως ATLAS (Automatically Tuned LA Software) βελτιστοποιημένη για την εκάστοτε αρχιτεκτονική.

Κάνοντας ένα βήμα πίσω όμως, θα δούμε ότι πολλές φορές χρειάζεται πολλαπλασιασμός διανύσματος με διάνυσμα, διανύσματος με πίνακα, και πίνακα με πίνακα. Για τις πράξεις αυτές υπάρχουν επίσης βελτιστοποιημένες βιβλιοθηκες, η BLAS (Basic Linear Algebra Subprograms), με loop unrolling που χωράει στην cache κλπ. Οποτε αν χρειάζεστε τις ανωτέρω πράξεις χρησιμοποιείστε τη BLAS.
http://www.netlib.org/blas/
Μια σημείωση μόνο. Αν χρειάζεστε να κάνετε πολλούς πολλαπλασιασμούς πινακών επί πινάκων---είναι από τους πιο αργούς αλγοριθμους που υπάρχει---ξανασκεφτείτε το σχεδιασμό του προγράμματός σας.

Φυσικά για πάρα πολλές αριθμητικές μεθόδους υπάρχουν ελεύθερες βιβλιοθήκες. Μερικά παραδείγματα:
  • Μη γραμμικές εξισώσεις και προβλήματα ελαχίστων τετραγώνων είναι η minpack: http://www.netlib.org/minpack/
  • Για επίλυση κανονικών διαφορικών εξισώσεων μια αρκετά στιβαρή βιβλιοθήκη είναι η odepack: http://www.netlib.org/odepack/ Χρησιμοποιείται και από το octave.
  • Μετασχηματισμός Fourier: http://www.netlib.org/fftpack/
  • Καμπύλες Bezier, splines γραμμένη από τον ίδιο τον de Boor: http://www.netlib.org/pppack/
Όπως θα παρατηρήσατε οι παραπάνω κώδικες φιλοξενούνται στην ιστοσελίδα http://www.netlib.org η οποία είναι μια τεράστια βάση. Απόλαυση να ψάχνεις.

Φυσικά δεν πρέπει να ξεχνάμε τη GNU βιβλιοθήκη http://www.gnu.org/software/gsl/ η οποία κάνει (σχεδόν) τα πάντα και συμφέρει (μιας και είναι δωρεάν)

Παιρνόντας στις μερικές διαφορικές εξισώσεις το τοπίο δυσκολεύει. Καταρχήν χρειάζεται να δημιουργηθεί το πλέγμα. Γι'αυτή τη δουλειά υπάρχουν τα gmsh και netgen. Μετά χρειάζονται οι βιβλιοθήκες για πεπερασμένους όγκους ή πεπερασμένα στοιχεία. Ενδεικτικά μερικές είναι libmesh, getfem, dealii και φυσικά άλλα πολλά.

Βαρετά όλα αυτά; Ας δούμε κάτι διαφορετικό. Αλγοριθμική παραγώγιση. Αντί να υπολογίσουμε τις παραγώγους με το χέρι ή με πεπερασμένες διαφορές (μην το κάνετε αυτό ποτέ, το πρόβλημα είναι ill-posed) μπορούμε να χρησιμοποιήσουμε αλγοριθμική ή αυτόματη παραγώγιση. Έχουμε ένα πρόγραμμα δηλαδή το οποίο δέχεται ως δεδομένα εισόδου τον κώδικα και βγάζει ως έξοδο ένα νέο κώδικα ο οποίος είναι η παράγωγος του αρχικού ως προς κάποιες μεταβλητές. Για παράδειγμα δείτε τα ADIFOR & ADOLC για Fortran & C αντίστοιχα, όπως και το Tapenade το οποίο είναι ελεύθερη υπηρεσία δικτύου.

Για γραφικά εκτός από OpenGL και τους canvaδες του εκάστοτε GUI, υπάρχουν οι βιβλιοθήκες PLPlot για διαγράμματα και άλλα συναφή και η VTK για 3D γραφικά πιο υψηλού επιπέδου από OpenGL.

H τακτική του mixed programming, δηλαδή να μπλέκονται πολλές γλώσσες προγραμματισμού στο ίδιο πρόγραμμα, δεν είναι άσχημη όταν υπάρχουν έτοιμες οι βιβλιοθήκες και απλώς γράφεται ένας wrapper. Για παράδειγμα, οι περισσότερες βιβλιοθήκες του numpy & scipy είναι wrappers για τις παραπάνω βιβλιοθήκες που ανέφερα.

Ακόμη πιο εξειδικευμένες βιβλιοθήκες είναι η http://gmplib.org/ για αριθμητική απεριόριστης ακρίβειας.

Αυτά προς το παρόν. Οποίος ενδιαφέρεται για κάτι πιο συγκεκριμένο ας το δηλώσει και κάτι μπορεί να γίνει. :geek:

Creative Commons License
Η εργασία υπάγεται στην άδεια Creative Commons Αναφορά-Παρόμοια διανομή 3.0 Ελλάδα

Re: Αριθμητικές μέθοδοι

ΔημοσίευσηΔημοσιεύτηκε: 24 Απρ 2010, 14:58
από gravity
"""""Παιρνόντας στις μερικές διαφορικές εξισώσεις το τοπίο δυσκολεύει. Καταρχήν χρειάζεται να δημιουργηθεί το πλέγμα. Γι'αυτή τη δουλειά υπάρχουν τα gmsh και netgen""" Tι κάνουν ακριβώς τα παραπάνω δυο προγράμματα?Βγάζουν πλέγμα για οποιοδήποτε γεωμετρία?

Re: Αριθμητικές μέθοδοι

ΔημοσίευσηΔημοσιεύτηκε: 24 Απρ 2010, 15:11
από Dimitris
Ακριβώς. Φυσικά το "οποιοδηποτε" σημαίνει ότι πρέπει να την εισάγεις εσύ στις μορφές αρχείων που υποστηρίζουν ή να τη σχεδιάσεις μέσα από το πρόγραμμα το ίδιο. Μετά το πλέγμα εξάγεται σε κάποια μορφή αρχείου για πλέγματα

Re: Αριθμητικές μέθοδοι

ΔημοσίευσηΔημοσιεύτηκε: 24 Απρ 2010, 15:28
από gravity
Ενδιαφέρον φαίνεται. Εγώ κάνω πεπερασμένες διαφορές και χρειάζομαι πλέγμα. Οπότε η ερώτηση τι είδους πλέγμα κάνουν. Φυσικά θα ψάξω και εγώ απλά ρωτάω γιατί πιθανό να γνωρίζεις. Εγώ συγκεκριμένα θέλω κυβικό πλέγμα και τώρα μάλλον σκέφτομαι να το παίρνω από το COMSOL αφού πρώτα το εγκαταστήσω sτα Ubuntu ή καλύτερα αφoύ προσπαθήσω να το εγκαταστήσω στο ubuntu.

Re: Αριθμητικές μέθοδοι

ΔημοσίευσηΔημοσιεύτηκε: 05 Δεκ 2010, 16:33
από Dimitris
Metanumerics

Στο παρελθόν είχα αναφερθεί σε κάποιο άλλο νήμα για μεταπρογραμματισμό (metaprogramming) και στο παρόν νήμα μιας και το θέμα είναι αριθμητικές μέθοδοι θα ήθελα να αναφερθώ στον όρο μετααριθμητική (metanumerics). Metanumerics είναι ένας σχετικά νέος κλάδος των εφαρμοσμένων μαθηματικών και της επιστημής των υπολογιστών που ασχολείται, κατ'αναλογία με το μεταπρογραμματισμό, με κώδικα που γράφει κώδικα. Στην περίπτωση αυτή βέβαια έχει να κάνει με κώδικα αριθμητικών μεθόδων. Ένα παράδειγμα metanumerics είχα επίσης παρουσιάσει σε ένα τεύχος του περιοδικού ubuntistas, τεύχος 9, για την αλγοριθμική ή αυτόματη παραγώγιση (algorithmic or automatic differentiation) και το πακέτο ADOL-C. Εδώ θα ήθελα να παρουσιάσω με περισσότερη λεπτομέρεια το θέμα metanumerics.

Ως metanumerics θα μπορούσαμε να ορίσουμε τον κώδικα που γράφει κώδικα αριθμητικών μεθόδων. Λίγο ασαφής ορισμός γι'αυτό ας δούμε ένα παράδειγμα. Έχουμε έναν κώδικα γραμμένο στη γλώσσα maxima το οποίο υπολογίζει το ολοκλήρωμα μιας συνάρτησης.
Κώδικας: Επιλογή όλων
integrate(x/(x^3 + 1), x);

Τίποτε δύσκολο. Ας υποθέσουμε ότι αυτό θέλουμε να το ενσωματώσουμε στον κώδικα C με τον οποίο δουλεύουμε. Η πιο απλή λύση είναι να προγραμματίσουμε τη λύση της ολοκλήρωσης σε μια συνάρτηση C. Το πρόβλημα που δημιουργείται είναι αν αλλάξουμε τη συνάρτηση ολοκλήρωσης στον κώδικα maxima, θα πρέπει να επαναπρογραμματίσουμε τη συνάρτηση C. Φυσικά στο απλό αυτό πρόβλημα θα μπορούσαμε να υλοποιήσουμε εναν αλγόριθμο αριθμητικής ολοκληρωσης κατ'ευθειαν στη C. Αλλα στη γενική περίπτωση θα μπορούσαμε να γράψουμε ένα script (θα παραμείνω στη χρήση του script γιατί είναι πιο εύκολο να γίνει σε μια interpreted γλώσσα) το οποίο θα παίρνει το αποτέλεσμα του maxima και θα δημιουργεί αυτόματα τον κώδικα C που χρειαζόμαστε, θα τον κάνει compile και έπειτα link στο C πρόγραμμά μας.

Το ADOL-C δέχεται ως δεδομένα εισόδου υπάρχωντα κώδικα τον οποίο τον επεξεργάζεται και δημιουργεί νέο κώδικα ο οποίος υπολογίζει την παράγωγο μιας μεταβλητής ως προς κάποια άλλη. Έπειτα ο νέος κώδικας μεταγλωττίζεται και συνδέεται στο πρόγραμμα όπου χρειάζεται η παράγωγος. Συνήθως είναι κάποιος κώδικας βελτιστοποίησης.

Ένα άλλο παράδειγμα προγράμματος που χρησιμοποιεί αυτή την τεχνική είναι το FEniCS. Το FEniCS είναι ένα υπερπακέτο, για python ή C++, το οποίο μπορεί να λύσει εξισώσεις στη variational form (βλ. calculus of variations = λογισμός των μεταβολών, μου λείπει η μετάφραση). Το FEniCS δημιουργεί κατά το runtime C++ κώδικα, ειδικό για το συγκεκριμένο πρόβλημα, ο οποίος μετά μεταγλωττίζεται, δημιουργείται η διεπιφάνεια για την python και τέλος εκτελείται.

Στο σημείο αυτό ίσως κάποιος πει ότι δημιουργείται ένα overhead και κάτι τέτοιο αληθεύει. Αλλά σε μεγάλα προβλήματα η γενικότητα του κώδικα δε συγκρίνεται με τη μικρή μείωση ταχύτητας. Μη ξεχνάμε ότι το κόστος του προγραμματιστή είναι πολύ μεγαλύτερο από το κόστος λειτουργίας του υπολογιστή. Να σημειώσω ότι το FEniCS λειτουργεί και παράλληλα με τη χρήση βιβλιοθηκών MPI.

Περνώντας πάλι στη γενική διατύπωση του ορισμού, μπορούμε να δώσουμε την εξής λύση σε πολλά προβλήματα. Ας υποθέσουμε ότι υπάρχει C/FORTRAN legacy code ο οποίος λύνει μια συγκεκριμένη κατηγορία προβλημάτων. Θα μπορούσε κάποιος να χρησιμοποιήσει μια interpreted γλώσσα προγραμματισμού για να δημιουργήσει case specific κώδικα C/FORTRAN ο οποίος θα μεταγλωττιστεί και θα ενσωματωθεί στο υπάρχων πρόγραμμα. Αυτό κάνουν οι βιβλιοθήκες boost στην python, οι lex/yacc για parsing, και πολλά άλλα εργαλεία.

Φυσικά ο προγραμματισμός αυτών των εργαλείων ίσως να απαιτήσει το χειρισμό strings (το οποίο για μένα είναι το χειρότερό μου) αλλά το τελικό αποτέλεσμα αξίζει τον κόπο.

Αλλάζοντας ελαφρώς θέμα θα ήθελα να αναφερθώ στις domain specific languages (DSL). Σε αντίθεση με τις γενικής χρήσης γλώσσες προγραμματισμού είναι γλώσσες που επιλύουν ένα συγκεκριμένο πρόβλημα και δεν ενδείκνυνται για προγραμματισμό όλων των τύπων προγραμμάτων. Για παράδειγμα το make αυτοματοποιεί τη μεταγλώττιση και σύνδεση των αρχείων ενός προγράμματος, το configure αυτοματοποιεί τη δημιουργεία των Makefiles. Συνήθως στόχος των DSL είναι η δημιουργία μιας γλώσσας με απλό συντακτικό.

Θα μπορούσαμε για παράδειγμα να γράψουμε ένα εργαλείο που θα έχει την εξής απλή σύνταξη:
Κώδικας: Επιλογή όλων
program {
  a -> b -> c;
  a: d -> e;
  b: f -> g;
  c: e -> g;
  f: d-> d;
}

το οποίο θα δημιουργεί τον εξής κώδικα C:
Κώδικας: Επιλογή όλων
/* main.c */
#include "file1.h"
#include "file2.h"
int main(){
a();
b();
c();
}

/* file1.c */
void a(){
d();
e();
}

void b(){
f();
g();
}

void c(){
e();
g();
}

/* file2.h */
void d(){
}

void e(){
}

void f(){
d();
d();
}

void g(){
}


και τα αντίστοιχα header files.

Αυτά για σήμερα. Καλό brainstorming και καλή διασκέδαση!!!