6

असीमित सीमा पर एक-आयामी एकीकरण के लिए संख्यात्मक एकीकरण (क्या संख्यात्मक विधि, और किस चाल का उपयोग करना है), जहां एकीकृत में एक या अधिक फ़ंक्शन 1d quantum harmonic oscillator तरंग कार्य हैं। अन्य लोगों के अलावा मैं हार्मोनिक दोलक आधार में कुछ समारोह के मैट्रिक्स तत्वों गणना करना चाहते हैं:क्वांटम हार्मोनिक ऑसीलेटर वेवफंक्शन के साथ संख्यात्मक एकीकरण कैसे करें?

फ़ाई n (x) = एन n एच n (एक्स) exp (-x/2)
जहां एच n (एक्स) Hermite polynomial

वी मीटर है, एन = \ पूर्णांक _ {- अनंत}^{} अनंत phi 012,मीटर (एक्स) वी (एक्स) phi n (x) dx

मामले में जहां अलग चौड़ाई के साथ क्वांटम हार्मोनिक wavefunctions देखते हैं में

इसके अलावा।

समस्या

, wavefunctions phi कि n (एक्स) oscillatory व्यवहार, जो बड़े n के लिए एक समस्या है, और GSL (जीएनयू वैज्ञानिक पुस्तकालय) से अनुकूली गॉस-Kronrod क्षेत्रकलन तरह एल्गोरिथ्म है की गणना करने के देर है और बड़ी त्रुटियां हैं।

+1

क्या यह एसओ पर कभी भी सबसे कठिन सवाल है? – MrTelly

+4

नहीं, यह केवल उस डोमेन से संबंधित है जो सबसे अधिक अपरिचित है, गूढ़! = कठिन। – Saem

+0

गॉस-लागुरेरे पहले ही [जीएसएल 2.3] (https://www.gnu.org/software/gsl/doc/html/integration.html) में पेश किए गए हैं – zmwang

उत्तर

8

एक अपूर्ण जवाब, क्योंकि मैं इस समय समय पर थोड़ा सा छोटा हूं; अगर अन्य तस्वीर को पूरा नहीं कर सकते हैं, तो मैं बाद में अधिक जानकारी प्रदान कर सकता हूं।

  1. जब भी और जहां भी संभव हो, तरंगों की ऑर्थोगोनैलिटी लागू करें। यह गणना की मात्रा में काफी कटौती करनी चाहिए।

  2. विश्लेषणात्मक रूप से जो कुछ भी आप कर सकते हैं। लिफ्ट स्थिरांक, भागों द्वारा एकीकृत विभाजित, जो भी हो। ब्याज के क्षेत्र को अलग करें; अधिकांश तरंगों बैंड-सीमित हैं, और ब्याज के क्षेत्र को कम करने से काम बचाने के लिए बहुत कुछ किया जाएगा।

  3. चतुर्भुज के लिए, आप शायद तरंगों को तीन टुकड़ों में विभाजित करना चाहते हैं और प्रत्येक को अलग से एकीकृत करना चाहते हैं: केंद्र में ऑसीलेटर बिट और दोनों तरफ घातीय-क्षय वाली पूंछ। यदि तरंगों का अजीब अजीब है, तो आप भाग्यशाली हो जाते हैं और पूंछ एक दूसरे को रद्द कर देंगे, जिसका अर्थ है कि आपको केवल केंद्र के बारे में चिंता करनी होगी। यहां तक ​​कि तरंगों के लिए, आपको केवल एक को एकीकृत करना होगा और इसे दोहराएं (समरूपता के लिए हुरे!)। अन्यथा, उच्च आदेश गॉस-लागुरेरे चौकोर नियम का उपयोग करके पूंछ को एकीकृत करें। आपको नियमों की गणना करनी पड़ सकती है; मुझे नहीं पता कि टेबल अच्छी गॉस-लागुरेरे नियमों की सूची में हैं, क्योंकि उनका उपयोग अक्सर नहीं किया जाता है। आप शायद त्रुटि व्यवहार को भी देखना चाहते हैं क्योंकि नियम में नोड्स की संख्या बढ़ जाती है; यह लंबे समय से रहा है क्योंकि मैंने गॉस-लागुरेरे नियमों का उपयोग किया था और मुझे याद नहीं है कि वे रनगे की घटना का प्रदर्शन करते हैं या नहीं। जो कुछ भी आपको पसंद है उसका उपयोग करके केंद्र भाग को एकीकृत करें; गॉस-क्रोन्रोड निश्चित रूप से एक ठोस विकल्प है, लेकिन फेजर क्वाड्रैचर भी होता है (जो कभी-कभी नोड्स की उच्च संख्या के लिए बेहतर होता है, जो एक ऑसीलेटर इंटीग्रैंड पर अच्छा काम कर सकता है) और यहां तक ​​कि ट्राइपोज़ाइडल नियम (जो कुछ oscillatory कार्यों के साथ आश्चर्यजनक सटीकता प्रदर्शित करता है))। एक उठाओ और इसे आजमाएं; यदि परिणाम खराब हैं, तो एक और तरीका एक शॉट दें।

SO पर कभी भी सबसे कठिन प्रश्न? मुश्किल से :)

0

मैं अभी इस बारे में किसी भी व्याख्या या योग्यता के लिए नहीं जा रहा हूं। यह कोड लिखा गया है और शायद गलत है। मुझे यह भी यकीन नहीं है कि यह वह कोड है जिसे मैं ढूंढ रहा था, मुझे बस याद है कि साल पहले मैंने यह समस्या की थी और मेरे अभिलेखागार को खोजने पर मुझे यह मिला। आपको खुद को आउटपुट प्लॉट करने की आवश्यकता होगी, कुछ निर्देश प्रदान किए जाते हैं। मैं कहूंगा कि अनंत सीमा पर एकीकरण एक समस्या है जिसे मैंने संबोधित किया और कोड के निष्पादन पर यह 'अनन्तता' पर गोल ऑफ त्रुटि को इंगित करता है (जो संख्यात्मक रूप से केवल बड़ा मतलब है)।

// compile g++ base.cc -lm 
#include <iostream> 
#include <cstdlib> 
#include <fstream> 
#include <math.h> 

using namespace std; 

int main() 
     { 
     double xmax,dfx,dx,x,hbar,k,dE,E,E_0,m,psi_0,psi_1,psi_2; 
     double w,num; 
     int n,temp,parity,order; 
     double last; 
     double propogator(double E,int parity); 
     double eigen(double E,int parity); 
     double f(double x, double psi, double dpsi); 
     double g(double x, double psi, double dpsi); 
     double rk4(double x, double psi, double dpsi, double E); 

     ofstream datas ("test.dat"); 

     E_0= 1.602189*pow(10.0,-19.0);// ev joules conversion 
     dE=E_0*.001; 
//w^2=k/m     v=1/2 k x^2    V=??? = E_0/xmax x^2  k--> 
//w=sqrt((2*E_0)/(m*xmax)); 
//E=(0+.5)*hbar*w; 

     cout << "Enter what energy level your looking for, as an (0,1,2...) INTEGER: "; 
     cin >> order; 

     E=0; 
     for (n=0; n<=order; n++) 
       { 
       parity=0; 
//if its even parity is 1 (true) 
       temp=n; 
       if ((n%2)==0) {parity=1; } 
       cout << "Energy " << n << " has these parameters: "; 
       E=eigen(E,parity); 
       if (n==order) 
         { 
         propogator(E,parity); 
         cout <<" The postive values of the wave function were written to sho.dat \n"; 
         cout <<" In order to plot the data should be reflected about the y-axis \n"; 
         cout <<" evenly for even energy levels and oddly for odd energy levels\n"; 
         } 
       E=E+dE; 
       } 
     } 

double propogator(double E,int parity) 
     { 
     ofstream datas ("sho.dat") ; 

     double hbar =1.054*pow(10.0,-34.0); 
     double m =9.109534*pow(10.0,-31.0); 
     double E_0= 1.602189*pow(10.0,-19.0); 
     double dx =pow(10.0,-10); 
     double xmax= 100*pow(10.0,-10.0)+dx; 
     double dE=E_0*.001; 
     double last=1; 
     double x=dx; 
     double psi_2=0.0; 
     double psi_0=0.0; 
     double psi_1=1.0; 
//  cout <<parity << " parity passsed \n"; 
     psi_0=0.0; 
     psi_1=1.0; 
     if (parity==1) 
       { 
       psi_0=1.0; 
       psi_1=m*(1.0/(hbar*hbar))* dx*dx*(0-E)+1 ; 
       } 

     do 
       { 
       datas << x << "\t" << psi_0 << "\n"; 
       psi_2=(2.0*m*(dx/hbar)*(dx/hbar)*(E_0*(x/xmax)*(x/xmax)-E)+2.0)*psi_1-psi_0; 
//cout << psi_1 << "=psi_1\n"; 
       psi_0=psi_1; 
       psi_1=psi_2; 
       x=x+dx; 
       } while (x<= xmax); 
//I return 666 as a dummy value sometimes to check the function has run 
     return 666; 
     } 


    double eigen(double E,int parity) 
     { 
     double hbar =1.054*pow(10.0,-34.0); 
     double m =9.109534*pow(10.0,-31.0); 
     double E_0= 1.602189*pow(10.0,-19.0); 
     double dx =pow(10.0,-10); 
     double xmax= 100*pow(10.0,-10.0)+dx; 
     double dE=E_0*.001; 
     double last=1; 
     double x=dx; 
     double psi_2=0.0; 
     double psi_0=0.0; 
     double psi_1=1.0; 
     do 
       { 
       psi_0=0.0; 
       psi_1=1.0; 

       if (parity==1) 
         {double psi_0=1.0; double psi_1=m*(1.0/(hbar*hbar))* dx*dx*(0-E)+1 ;} 
       x=dx; 
       do 
         { 
         psi_2=(2.0*m*(dx/hbar)*(dx/hbar)*(E_0*(x/xmax)*(x/xmax)-E)+2.0)*psi_1-psi_0; 
         psi_0=psi_1; 
         psi_1=psi_2; 
         x=x+dx; 
         } while (x<= xmax); 


       if (sqrt(psi_2*psi_2)<=1.0*pow(10.0,-3.0)) 
         { 
         cout << E << " is an eigen energy and " << psi_2 << " is psi of 'infinity' \n"; 
         return E; 
         } 
       else 
         { 
         if ((last >0.0 && psi_2<0.0) ||(psi_2>0.0 && last<0.0)) 
           { 
           E=E-dE; 
           dE=dE/10.0; 
           } 
         } 
       last=psi_2; 
       E=E+dE; 
       } while (E<=E_0); 
     } 

यदि यह कोड सही, गलत, रोचक लगता है या आपके पास विशिष्ट प्रश्न पूछते हैं और मैं उनका उत्तर दूंगा।

4

मैं कुछ अन्य बातें की सलाह देते हैं:

  1. एक परिमित डोमेन पर समारोह बदलने एकीकरण अधिक प्रबंधनीय बनाने की कोशिश करें।
  2. जहां संभव हो समरूपता का उपयोग करें - इसे नकारात्मक अनंतता से शून्य और शून्य तक दो इंटीग्रल के योग में विभाजित करें और देखें कि फ़ंक्शन समरूपता या विरोधी सममित है या नहीं। यह आपके कंप्यूटिंग को आसान बना सकता है।
  3. Gauss-Laguerre quadrature पर देखें और देखें कि यह आपकी मदद कर सकता है या नहीं।
0

मैं भौतिकी में एक छात्र प्रमुख हूं, और मुझे भी समस्या का सामना करना पड़ा। इन दिनों मैं इस सवाल के बारे में सोचता रहता हूं और अपना उत्तर प्राप्त करता हूं। मुझे लगता है कि यह आपको इस प्रश्न को हल करने में मदद कर सकता है।

1. जीएसएल में, फ़ंक्शन हैं जो आपको ऑसीलेटर फ़ंक्शन को एकीकृत करने में मदद कर सकते हैं - qawo & qawf। शायद आप एक मान सेट कर सकते हैं, । और एकीकरण को टॉव भागों में विभाजित किया जा सकता है, [0, ] और [, pos_infinity]। पहले अंतराल में, आप किसी भी जीएसएल एकीकरण फ़ंक्शन का उपयोग कर सकते हैं, और दूसरे अंतराल में, आप qawo या qawf का उपयोग कर सकते हैं।

2.Or आप एक ऊपरी सीमा के समारोह को एकीकृत कर सकते, , कि [0, ] में एकीकृत है। इसलिए एकीकरण को गॉस लीजेंड्री विधि का उपयोग करके गणना की जा सकती है, और यह जीएसएल में प्रदान की जाती है। यद्यपि वास्तविक मूल्य और गणना मूल्य के बीच कुछ अंतर हो सकता है, लेकिन यदि आप बी ठीक से सेट करते हैं, तो अंतर को उपेक्षित किया जा सकता है। जब तक आप चाहते हैं कि सटीकता से अंतर कम है। और जीएसएल फ़ंक्शन का उपयोग करने वाला यह तरीका केवल एक बार बुलाया जाता है और कई बार उपयोग कर सकता है, क्योंकि रिटर्न वैल्यू पॉइंट और उसके इसी वजन का होता है, और एकीकरण केवल एफ (xi) * wi का योग होता है, अधिक जानकारी के लिए आप गॉस लीजेंड्रे को खोज सकते हैं विकिपीडिया पर चौकोर। एकाधिक और अतिरिक्त ऑपरेशन एकीकरण से बहुत तेज है।

3. एक ऐसा कार्य भी है जो अनंत क्षेत्र एकीकरण की गणना कर सकता है - qagi, आप इसे जीएसएल-उपयोगकर्ता की मार्गदर्शिका में खोज सकते हैं। लेकिन इसे हर बार बुलाया जाता है जब आपको एकीकरण की गणना करने की आवश्यकता होती है, और इससे कुछ समय लग सकता है, लेकिन मुझे यकीन नहीं है कि यह आपके कार्यक्रम में कितनी देर तक उपयोग करेगा।

मैं सुझाव देता हूं कि मैंने प्रस्तावित नंबर 2 पसंद किया है।

0

आप n से भी कम समय = 100 हार्मोनिक दोलक कार्यों के साथ काम करने जा रहे हैं, तो आप कोशिश करना चाहते हो सकता है:

http://www.mymathlib.com/quadrature/gauss_hermite.html

कार्यक्रम 100 शून्य और वजन (साथ गॉस-हर्मिट क्षेत्रकलन के माध्यम से एक अभिन्न गणना करता है H_100 के शून्य)। एक बार जब आप Hermite_100 पर जाते हैं तो इंटीग्रल सटीक नहीं होते हैं।

इस एकीकरण विधि का उपयोग करके मैंने एक प्रोग्राम लिखा है जो आप गणना करना चाहते हैं और यह काफी अच्छी तरह से काम करता है। इसके अलावा, Hermite-polynomial शून्य के एसिम्प्टोटिक रूप का उपयोग कर एन = 100 से आगे जाने का एक तरीका हो सकता है लेकिन मैंने इसे नहीं देखा है।

संबंधित मुद्दे