(*==========================================================================*) (* *) (* PASCAL-XSC - PROGRAMM PIVP9906 990531 *) (* *) (* PROGRAMM ZUR L™SUNG DES ANFANGSWERTPROBLEMS ZUR DGL *) (* *) (* y^(n) = \sum_(i=0)^(n-1) p_i(x) y^(i) + p(x), p(x), p_i(x) Polynome *) (* *) (* POTENZREIHENANSATZ: y(x) = \sum_0^infty a_k x^k *) (* L™SUNG MIT VARIABLER TAYLORENTWICKLUNG UND RESTGLIEDABSCHŽTZUNG *) (* RECHNUNG MIT ADAPTIVER GENAUIGKEIT *) (* *) (* HIER : HAUPTPROGRAMM *) (* EINLESEN DES ANFANGSWERTPROBLEMS *) (* ZERLEGEN DES DEFINITIONSINTERVALLS IN GEEIGNETE ABSCHNITTE *) (* AUFRUF DER PROZEDUR IVP_POLY IM MODUL PIVP_MPS *) (* *) (*==========================================================================*) PROGRAM PIVP9906( INPUT , OUTPUT ); USE iostd , x_real, i_ari , mv_ari , mvi_ari , lss , ilss , service , mps_arip , pivp_mps , timer; CONST MMAX = 4; (* MAXIMALER GRAD DER POLYNOME *) KMAX = 2000; (* MAXIMALE ORDNUNG DER TAYLORENTWICKLUNG VON Y *) (* T_IND : TEST-INDEX FšR ABBRUCHKRITERIUM IM UNTERPROGRAMM IVP_POLY *) (* 0 FšR FUNKTION SELBST, J FšR J-TE ABLEITUNG *) VAR N , P , SCHRITTE , T_IND , VARIANTEN : INTEGER; ERR , INDEKS , SUBINTERVALS , INHDGL : INTEGER; INHOM : BOOLEAN; E_ABS , CPU_TIME : REAL; AKT_DATUM , AWP_STRING : STRING; LOESUNG : STRING; DATEIN , OUTFILE , CONSOLE : TEXT; (*--------------------------------------------------------------------------*) PROCEDURE EINLESEN_MAIN; BEGIN RESET( INPUT , 'awp.ein' ); REWRITE( OUTPUT , 'awp.aus' ); REWRITE( OUTFILE , 'outfile.awp' ); REWRITE( CONSOLE , 'con' ); (* SYSTEMDATUM EINLESEN AUS FILE *) RESET( DATEIN , 'aktdatum.txt' ); READLN( DATEIN , AKT_DATUM ); (* FUNKTIONSNAME EINLESEN *) READLN; READLN; READLN( AWP_STRING ); READLN; READLN( LOESUNG ); READLN; (* ORDNUNG DER DGL EINLESEN *) READLN( N ); (* HOMOGENE ODER INHOMOGENE DGL? *) READLN( INHDGL ); IF INHDGL = 0 THEN INHOM := FALSE ELSE INHOM := TRUE; END; (*--------------------------------------------------------------------------*) PROCEDURE IVP_HAUPT( N : INTEGER ); (* VARIABLEN : *) (* TRB : POLYNOMKOEFFIZIENTEN *) (* X_END : RECHTE GRENZE DES DEFINITIONSINTERVALLS DER DGL *) VAR I , J , K : INTEGER; M_AKT_MAX , AW_REAL , INH_GRAD : INTEGER; X_AKT , X_END : REAL; VRM : INTEGERVECTOR[0..N]; VRY , VRNULL : RVECTOR[0..N-1]; VR_PART1 , VR_PART4 , VR_PART5 : RVECTOR[0..N-1]; VIY , VIZ , VI_INH : IVECTOR[0..N-1]; VI_PART1 , VI_PART4 , VI_PART5 : IVECTOR[0..N-1]; VINULL , VIY1 : IVECTOR[0..N-1]; VIY2 , VIY3 , VIY4 , VIY5 : IVECTOR[0..N-1]; VIU2 , VIU3 , VIU4 , VIU5 : IVECTOR[0..N-1]; MRID : RMATRIX[0..N-1,0..N-1]; MRB , MRB2 , MRB3 , MRB4 , MRB5 : RMATRIX[0..N-1,0..N-1]; TRB : RMATRIX[0..N,0..MMAX]; MIA , MIBINV2 , MIBINV3 : IMATRIX[0..N-1,0..N-1]; MIBINV4 , MIBINV5 : IMATRIX[0..N-1,0..N-1]; FAK : RVECTOR[0..N]; PRODINDEX : INTEGER; PROD : RMATRIX[0..KMAX,0..N]; IPROD : IMATRIX[0..KMAX,0..N]; (*--------------------------------------------------------------------------*) PROCEDURE EINLESEN; VAR I , J : INTEGER; BEGIN (* ORDNUNG DER POLYNOME SOWIE POLYNOMKOEFFIZIENTEN EINLESEN *) FOR I := N - 1 DOWNTO 0 DO BEGIN READLN( VRM[I] ); FOR J := VRM[I] DOWNTO 0 DO READLN( TRB[I,J] ); READLN; END; (* ORDNUNG UND POLYNOMKOEFFIZIENTEN DER INHOMOGENITŽT EINLESEN *) (* SPEICHERN MIT ERSTEM INDEX N, DER SONST UNBELEGT IST *) READLN( VRM[N] ); INH_GRAD := VRM[N]; FOR J := VRM[N] DOWNTO 0 DO READLN( TRB[N,J] ); READLN; (* ANFANGSWERTE EINLESEN *) READLN( AW_REAL ); IF AW_REAL = 0 THEN (* REELLE ANFANGSWERTE *) FOR I := 0 TO N - 1 DO BEGIN READLN( VRY[I] ); VIY[I] := VRY[I]; END ELSE FOR I := 0 TO N - 1 DO READLN( VIY[I] ); (* RECHTE GRENZE DES DEFINITIONSINTERVALLS DER DGL EINLESEN *) READLN( X_END ); (* ABSOLUTEN GENAUIGKEITSPARAMETER EINLESEN *) (* NUR FšR UNTERLAUFBEREICH INTERESSANT *) READLN( E_ABS ); (* TEST-INDEX EINLESEN *) (* T_IND-te ABLEITUNG WIRD AUF GENAUIGKEIT šBERPRšFT *) READLN( T_IND ); (* SCHRITTE EINLESEN (1 FšR EINEN SCHRITT, 0 FšR GESTEUERTE SCHRITTWEITE *) READLN( SCHRITTE ); (* TEST-ORDNUNG DER TAYLOR-ENTWICKLUNG EINLESEN *) (* FALLS DIE KONTRAKTIONSBEDINGUNG NACH P SUMMANDEN *) (* NICHT ERFšLLT IST, WIRD DAS INTERVALL HALBIERT *) (* NUR BEI GESTEUERTER SCHRITTWEITE BENUTZT *) READLN( P ); (* EINSCHLIESSUNGSVARIANTEN EINLESEN *) READLN( VARIANTEN ); END; (*--------------------------------------------------------------------------*) PROCEDURE WRITE_HEAD; VAR I : INTEGER; BEGIN WRITE( 'Verifizierte L”sung des Anfangswertproblems ' ); WRITELN( ' ' , AKT_DATUM ); WRITELN; WRITELN( AWP_STRING ); WRITELN; FOR I := 0 TO N - 1 DO WRITELN( 'y^(' , I : 1 , ') (0) = ' , VIY[I].SUP : 12 : 8 , ' ' ); WRITELN; WRITELN( 'im Intervall [ 0 , ' , X_END , ' ]' ); WRITELN; WRITELN( LOESUNG ); WRITELN; WRITELN( 'Programm PIVP9906 vom 31.05.99' ); WRITELN; WRITELN( 'Potenzreihenansatz mit adaptiver Genauigkeit.' ); IF SCHRITTE = 1 THEN WRITELN( 'Rechnung mit einem Integrationsschritt.' ) ELSE BEGIN WRITE( 'Schrittweitensteuerung mit Kontraktionsbedingung nach ' ); WRITELN( P : 1 , ' Koeffizienten.' ); END; WRITELN( 'Die Kontraktionsbedingung wird fr die (n-1)-te Ableitung geprft.' ); WRITELN; WRITELN( 'Abbruchkriterium: E_abs = ' , E_ABS , '.' ); WRITE( ' šberprft wurde die Genauigkeit der ' ); WRITELN( T_IND : 1 , '. Ableitung.' ); WRITELN; END; (*--------------------------------------------------------------------------*) PROCEDURE QR( VAR MRA : RMATRIX ); (* QR-ZERLEGUNG DER MATRIX MRA *) (* AUSGABE IM PARAMETER MRA : ORTHOGONALE MATRIX MRQ *) VAR I , J , K , NL , NU : INTEGER; MRQ : RMATRIX[LBOUND(MRA)..UBOUND(MRA),LBOUND(MRA)..UBOUND(MRA)]; B,C : RVECTOR[LBOUND(MRA)..UBOUND(MRA)]; S : REAL; BEGIN NL := LBOUND(MRA); NU := UBOUND(MRA); MRQ := ID( MRQ ); FOR K := NL TO NU-1 DO BEGIN S := 0; FOR I := K TO NU DO BEGIN B[I] := MRA[I,K]; IF S < ABS( B[I] ) THEN S := ABS( B[I] ); END; IF S > 0 THEN BEGIN S := SQRT( B * B ); IF MRA[K,K] >= 0 THEN S := -S; B[K] := MRA[K,K] - S; S := 1 / ( S * B[K] ); C := S * B; FOR J := K + 1 TO NU DO BEGIN S := C * RVECTOR( MRA[*,J] ); FOR I := K TO NU DO MRA[I,J] := MRA[I,J] + B[I]*S; END; FOR I := NL TO NU DO BEGIN S := RVECTOR( MRQ[I] ) * C; FOR J := K TO NU DO MRQ[I,J] := MRQ[I,J] + S * B[J]; END; END; B[K] := 0; C[K] := 0; END; MRA := MRQ; END; (*--------------------------------------------------------------------------*) PROCEDURE QR_ZERL( VAR PHI : RMATRIX; REST : IVECTOR ); (* VERANLASST QR-ZERLEGUNG VON PHI NACH SORTIEREN DER SPALTEN VON PHI *) (* GLOBALE VARIABLE : N *) VAR I , J , K : INTEGER; HV , LAENGE : RVECTOR[0..N-1]; HR : REAL; BEGIN FOR I := 0 TO N - 1 DO LAENGE[I] := SQR( DIAM( REST[I] ) ) * ( RVECTOR( PHI[*,I] ) * RVECTOR( PHI[*,I] ) ); FOR J := 0 TO N - 2 DO FOR K := J + 1 TO N - 1 DO IF LAENGE[J] < LAENGE[K] THEN BEGIN HR := LAENGE[K]; LAENGE[K] := LAENGE[J]; LAENGE[J] := HR; HV := RVECTOR(PHI[*,K]); PHI[*,K] := PHI[*,J]; PHI[*,J] := HV; END; QR( PHI ); (* PHI WIRD MIT Q AUS DER QR-ZERLEGUNG VON PHI šBERSCHRIEBEN *) END; (*--------------------------------------------------------------------------*) PROCEDURE FUNSYS( VRY: RVECTOR; VAR MRB : RMATRIX; VAR INDEKS : INTEGER ); (* LOKALES FUNDAMENTALSYSTEM *) VAR I : INTEGER; MAX : REAL; BEGIN MAX := ABS( VRY[0] ); INDEKS := 0; FOR I := 1 TO N - 1 DO IF ABS( VRY[I] ) > MAX THEN BEGIN MAX := ABS( VRY[I] ); INDEKS := I; END; MRB := ID( MRB ); MRB[*,INDEKS] := VRY; END; (*--------------------------------------------------------------------------*) FUNCTION FUNSYSLGS( VAR MRB: RMATRIX; INDEKS: INTEGER; VAR VIY: IVECTOR ) : IVECTOR[0..N-1]; (* LOEST LGS MIT SPEZIELLER MATRIX *) VAR I : INTEGER; VIZ : IVECTOR[0..N-1]; BEGIN FOR I := 0 TO N - 1 DO IF I = INDEKS THEN VIZ[I] := VIY[I] / MRB[INDEKS,INDEKS] ELSE VIZ[I] := VIY[I] - VIY[INDEKS] * MRB[I,INDEKS] / MRB[INDEKS,INDEKS]; FUNSYSLGS := VIZ; END; (*--------------------------------------------------------------------------*) PROCEDURE IVP_STEP( VAR X_AKT : REAL; H : REAL; VARIANTEN: INTEGER ); VAR I , J , K , M : INTEGER; MPX_AKT_POT : MPVECTOR[0..MMAX]; TIC : MPMATRIX[0..N,0..MMAX]; (* LOKALE POLYNOMKOEFFIZIENTEN *) MP_TRP : MPVECTOR[0..MMAX]; MPH : MPVECTOR[0..N+MMAX]; VRY : RVECTOR[0..N-1]; (*----------------*) FUNCTION DECREASE( K : INTEGER ) : BOOLEAN; (* SCHŽTZT AB, OB DIE KOEFFIZIENTEN AB DEM K-TEN BETRAGSMŽSSIG FALLEN *) (* DAS WACHSTUM DER KOEFFIZIENTEN WIRD GERINGFšGIG UNTERSCHŽTZT *) VAR I , J , L : INTEGER; KOEFSUM : REAL; LESSONE : BOOLEAN; TRP : RVECTOR[0..MMAX]; TRC : RMATRIX[0..N,0..MMAX]; (* LOKALE POLYNOMKOEFFIZIENTEN *) BEGIN LESSONE := FALSE; KOEFSUM := 0; IF X_AKT = 0 THEN FOR I := N - 1 DOWNTO 0 DO FOR J := 0 TO VRM[I] DO TRC[I,J] := TRB[I,J] ELSE BEGIN FOR I := 0 TO N DO IF VRM[I] >= 0 THEN BEGIN FOR J := 0 TO VRM[I] DO TRP[J] := TRB[I,J]; FOR J := 0 TO VRM[I] DO BEGIN TRC[I,J] := TRP[0]; FOR L := 1 TO VRM[I] - J DO TRC[I,J] := TRC[I,J] + ( X_AKT ** L ) * TRP[L]; TRC[I,J] := TRC[I,J] / FAC( J ); FOR L := 0 TO VRM[I] - J - 1 DO TRP[L] := TRP[L+1] * ( L + 1 ); END; END; END; FOR I := N - 1 DOWNTO 0 DO BEGIN FOR J := 0 TO VRM[I] DO KOEFSUM := KOEFSUM + ABS( TRC[I,J] ) * ( H ** J ); KOEFSUM := K * KOEFSUM / H; END; KOEFSUM := KOEFSUM * ( ( H / K ) ** ( N + 1 ) ); {WRITELN( OUTFILE ); WRITELN( OUTFILE, 'KOEFSUM nach ',K:1,' Koeffizienten ca. ',KOEFSUM:10:4 ); FLUSH( OUTFILE );} IF KOEFSUM < 1 THEN LESSONE := TRUE; DECREASE := LESSONE; END; (*----------------*) (*-------------------------------------------------------------------------*) BEGIN (* IVP_STEP *) (* ANFANGSWERTE VORBESETZEN *) FOR I := 0 TO N DO FOR J := 0 TO MMAX DO TIC[I,J].PREC := 1; M_AKT_MAX := VRM[0]; FOR I := 1 TO N DO IF VRM[I] > M_AKT_MAX THEN M_AKT_MAX := VRM[I]; (* POTENZEN VON H BEREITSTELLEN *) MPH[0].PREC := 1; MPH[0].STAG := 1; MPH[0].INT := 0; MPH[1].PREC := 1; MPH[1].STAG := H; MPH[1].INT := 0; FOR I := 2 TO N + MMAX DO MPH[I] := MPH[I-1] *> X_END; IF ( SCHRITTE = 1 ) OR ( DECREASE( P ) ) THEN (* NUR BEI REELLEN ANFANGSWERTEN IST EINSCHRITTVERFAHREN M™GLICH *) (* BEI INTERVALL-ANFANGSWERTEN WIRD FUNDAMENTALSYSTEM EINGESCHLOSSEN *) IF ( AW_REAL = 0 ) AND ( SCHRITTE = 1 ) THEN BEGIN FOR I := 0 TO N DO FOR J := 0 TO MMAX DO BEGIN TIC[I,J].PREC := 1; TIC[I,J].STAG := TRB[I,J]; TIC[I,J].INT := 0; END; (* POLYNOMKOEFFIZIENTEN DER INHOMOGENITŽT *) (* MIT POTENZEN VON MPH MULTIPLIZIEREN *) FOR J := VRM[N] DOWNTO 0 DO TIC[N,J] := TIC[N,J] *> MPH[N+J]; IVP_POLY( N , M_AKT_MAX , VRM , TIC , MID( VIY ) , H , E_ABS , T_IND , FAK , PRODINDEX , PROD , IPROD , VIZ , OUTFILE , CONSOLE ); VIY := VIZ; WRITELN; WRITELN( 'Einschlieáungen an der Stelle ' , H , ':' ); WRITELN; FOR I := 0 TO N - 1 DO BEGIN WRITE( 'y^(' , I : 1 , ') (' , H : 10 : 9 , ') î ' ); IWRITELN( VIY[I] ); END; END ELSE BEGIN IF X_AKT > 0 THEN BEGIN (* POLYNOME NACH ( X - X_AKT ) ENTWICKELN *) (* ZUERST POTENZEN VON X_AKT BEREITSTELLEN *) MPX_AKT_POT[0].PREC := 1; MPX_AKT_POT[0].STAG := 1; MPX_AKT_POT[0].INT := 0; MPX_AKT_POT[1].PREC := 1; MPX_AKT_POT[1].STAG := X_AKT; MPX_AKT_POT[1].INT := 0; FOR K := 2 TO M_AKT_MAX DO BEGIN MPX_AKT_POT[K] := MPX_AKT_POT[K-1] *> X_AKT; END; (* POLYNOME NACH ( X - X_AKT ) ENTWICKELN *) FOR I := 0 TO N DO IF VRM[I] >= 0 THEN BEGIN FOR J := 0 TO VRM[I] DO BEGIN MP_TRP[J].STAG := TRB[I,J]; MP_TRP[J].INT := 0; MP_TRP[J].PREC := 1; END; FOR J := 0 TO VRM[I] DO BEGIN TIC[I,J] := MP_TRP[0]; FOR K := 1 TO VRM[I] - J DO BEGIN TIC[I,J] := TIC[I,J] + ( MPX_AKT_POT[K] *< MP_TRP[K] ); END; TIC[I,J] := TIC[I,J] / FAC( J ); FOR K := 0 TO VRM[I] - J - 1 DO MP_TRP[K] := MP_TRP[K+1] *> ( K + 1 ); END; END; END ELSE (* POLYNOMKOEFFIZIENTEN UNVERŽNDERT IN TIC SPEICHERN *) FOR I := 0 TO N DO FOR J := 0 TO MMAX DO BEGIN TIC[I,J].PREC := 1; TIC[I,J].STAG := TRB[I,J]; TIC[I,J].INT := 0; END; (* POLYNOMKOEFFIZIENTEN DER INHOMOGENITŽT *) (* MIT POTENZEN VON MPH MULTIPLIZIEREN *) FOR J := VRM[N] DOWNTO 0 DO TIC[N,J] := TIC[N,J] *> MPH[N+J]; X_AKT := X_AKT + H; WRITELN( CONSOLE , 'X_AKT = ' , X_AKT ); FLUSH( CONSOLE ); WRITELN; WRITELN( 'Einschlieáungen an der Stelle ' , X_AKT , ':' ); (* AWP-L™SER FšR FUNDAMENTALSYSTEM AUFRUFEN *) (* VARIANTE 1: MITTELWERTVERFAHREN *) IF VARIANTEN > 15 THEN BEGIN VARIANTEN := VARIANTEN - 16; (* MITTELWERTL™SUNG *) IVP_POLY( N , M_AKT_MAX , VRM , TIC , VR_PART1 , H , E_ABS , T_IND , FAK , PRODINDEX , PROD , IPROD , VI_PART1 , OUTFILE , CONSOLE ); IF VIY1 <> VINULL THEN BEGIN (* HOMOGENE DGL; RECHTE SEITE DURCH VRM[N] AUSBLENDEN *) VRM[N] := -1; FOR J := 0 TO N - 1 DO BEGIN VRY := MRID[J]; IVP_POLY( N , M_AKT_MAX , VRM , TIC , VRY , H , E_ABS , T_IND , FAK , PRODINDEX , PROD , IPROD , VIZ , OUTFILE , CONSOLE ); MIA[*,J] := VIZ; END; VRM[N] := INH_GRAD; END; VIY1 := MIA * VIY1 + VI_PART1; WRITELN; WRITELN( 'Mittelwertverfahren:' ); WRITELN; FOR I := 0 TO N - 1 DO BEGIN WRITE( 'y^(' , I : 1 , ') (' , X_AKT : 10 : 9 , ') î ' ); IWRITELN( VIY1[I] ); END; VR_PART1 := MID( VIY1 ); VIY1 := VIY1 - VR_PART1; WRITELN; END; (* VARIANTE 2: PARALELLEPIPEDE *) (* MITTELPUNKTSMATRIX ALS FUNDAMENTALSYSTEM *) (* INHOMOGENE DGL MIT HOMOGENEN ABD *) (* FšR INHOMOGENE DGLN FRAGWšRDIG *) IF VARIANTEN > 7 THEN BEGIN VARIANTEN := VARIANTEN - 8; IF INHOM THEN (* PARTIKULŽRL™SUNG MIT HOMOGENEN ANFANGSWERTEN *) IVP_POLY( N , M_AKT_MAX , VRM , TIC , VRNULL , H , E_ABS , T_IND , FAK , PRODINDEX , PROD , IPROD , VI_INH , OUTFILE , CONSOLE ); (* HOMOGENE DGL; RECHTE SEITE DURCH VRM[N] AUSBLENDEN *) VRM[N] := -1; FOR J := 0 TO N - 1 DO BEGIN VRY := MRB2[*,J]; IVP_POLY( N , M_AKT_MAX , VRM , TIC , VRY , H , E_ABS , T_IND , FAK , PRODINDEX , PROD , IPROD , VIZ , OUTFILE , CONSOLE ); MIA[*,J] := VIZ; END; VRM[N] := INH_GRAD; (* PARALLEL-EPIPED-AUSWERTUNG *) MRB2 := MID( MIA ); INV( MRB2 , MIBINV2 , K ); VIU2 := ( MIBINV2 * MIA ) * VIU2 + MIBINV2 * VI_INH; VIY2 := MRB2 * VIU2; WRITELN; WRITELN( 'Fortsetzung mit Parallelepipeden:' ); WRITELN; FOR I := 0 TO N - 1 DO BEGIN WRITE( 'y^(' , I : 1 , ') (' , X_AKT : 10 : 9 , ') î ' ); IWRITELN( VIY2[I] ); END; WRITELN; END; (* VARIANTE 3: QR-ZERLEGUNG *) (* INHOMOGENE DGL MIT HOMOGENEN ABD *) (* FšR INHOMOGENE DGLN FRAGWšRDIG *) IF VARIANTEN > 3 THEN BEGIN VARIANTEN := VARIANTEN - 4; IF INHOM THEN (* PARTIKULŽRL™SUNG MIT HOMOGENEN ANFANGSWERTEN *) IVP_POLY( N , M_AKT_MAX , VRM , TIC , VRNULL , H , E_ABS , T_IND , FAK , PRODINDEX , PROD , IPROD , VI_INH , OUTFILE , CONSOLE ); (* HOMOGENE DGL; RECHTE SEITE DURCH VRM[N] AUSBLENDEN *) VRM[N] := -1; FOR J := 0 TO N - 1 DO BEGIN VRY := MRB3[*,J]; IVP_POLY( N , M_AKT_MAX , VRM , TIC , VRY , H , E_ABS , T_IND , FAK , PRODINDEX , PROD , IPROD , VIZ , OUTFILE , CONSOLE ); MIA[*,J] := VIZ; END; VRM[N] := INH_GRAD; (* QR-AUSWERTUNG *) MRB3 := MID( MIA ); QR_ZERL( MRB3 , VIU3 ); INV( MRB3 , MIBINV3 , K ); VIU3 := ( MIBINV3 * MIA ) * VIU3 + MIBINV3 * VI_INH; VIY3 := MRB3 * VIU3; WRITELN; WRITELN( 'Fortsetzung mit QR-Zerlegung:' ); WRITELN; FOR I := 0 TO N - 1 DO BEGIN WRITE( 'y^(' , I : 1 , ') (' , X_AKT : 10 : 9 , ') î ' ); IWRITELN( VIY3[I] ); END; WRITELN; END; (* VARIANTE 4: ERWEITERTES MITTELWERTVERFAHREN *) (* MITTELPUNKTSMATRIX ALS FUNDAMENTALSYSTEM *) IF VARIANTEN > 1 THEN BEGIN VARIANTEN := VARIANTEN - 2; (* MITTELWERTL™SUNG *) IVP_POLY( N , M_AKT_MAX , VRM , TIC , VR_PART4 , H , E_ABS , T_IND , FAK , PRODINDEX , PROD , IPROD , VI_PART4 , OUTFILE , CONSOLE ); IF VIY4 <> VINULL THEN BEGIN (* HOMOGENE DGL; RECHTE SEITE DURCH VRM[N] AUSBLENDEN *) VRM[N] := -1; FOR J := 0 TO N - 1 DO BEGIN VRY := MRB4[*,J]; IVP_POLY( N , M_AKT_MAX , VRM , TIC , VRY , H , E_ABS , T_IND , FAK , PRODINDEX , PROD , IPROD , VIZ , OUTFILE , CONSOLE ); MIA[*,J] := VIZ; END; VRM[N] := INH_GRAD; END; (* PARALLEL-EPIPED-AUSWERTUNG *) MRB4 := MID( MIA ); INV( MRB4 , MIBINV4 , K ); VIU4 := ( MIBINV4 * MIA ) * VIU4 + MIBINV4 * VI_PART4; VIY4 := MRB4 * VIU4; WRITELN; WRITELN( 'Erweitertes Mittelwertverfahren mit Mittelpunktsmatrix:' ); WRITELN; FOR I := 0 TO N - 1 DO BEGIN WRITE( 'y^(' , I : 1 , ') (' , X_AKT : 10 : 9 , ') î ' ); IWRITELN( VIY4[I] ); END; VR_PART4 := MID( VIY4 ); VIU4 := VIU4 - MIBINV4 * VR_PART4; WRITELN; END; (* VARIANTE 5: ERWEITERTES MITTELWERTVERFAHREN MIT QR-ZERLEGUNG *) IF VARIANTEN > 0 THEN BEGIN (* MITTELWERTL™SUNG *) IVP_POLY( N , M_AKT_MAX , VRM , TIC , VR_PART5 , H , E_ABS , T_IND , FAK , PRODINDEX , PROD , IPROD , VI_PART5 , OUTFILE , CONSOLE ); IF VIY5 <> VINULL THEN BEGIN (* HOMOGENE DGL; RECHTE SEITE DURCH VRM[N] AUSBLENDEN *) VRM[N] := -1; FOR J := 0 TO N - 1 DO BEGIN VRY := MRB5[*,J]; IVP_POLY( N , M_AKT_MAX , VRM , TIC , VRY , H , E_ABS , T_IND , FAK , PRODINDEX , PROD , IPROD , VIZ , OUTFILE , CONSOLE ); MIA[*,J] := VIZ; END; VRM[N] := INH_GRAD; END; (* QR-AUSWERTUNG *) MRB5 := MID( MIA ); QR_ZERL( MRB5 , VIU5 ); INV( MRB5 , MIBINV5 , K ); VIU5 := ( MIBINV5 * MIA ) * VIU5 + MIBINV5 * VI_PART5; VIY5 := MRB5 * VIU5; WRITELN; WRITELN( 'Erweitertes Mittelwertverfahren mit QR-Zerlegung:' ); WRITELN; FOR I := 0 TO N - 1 DO BEGIN WRITE( 'y^(' , I : 1 , ') (' , X_AKT : 10 : 9 , ') î ' ); IWRITELN( VIY5[I] ); END; VR_PART5 := MID( VIY5 ); VIU5 := VIU5 - MIBINV5 * VR_PART5; WRITELN; END; END ELSE BEGIN (* ANZAHL DER TEILINTERVALLE ERH™HEN *) SUBINTERVALS := SUBINTERVALS + 1; IVP_STEP( X_AKT , H / 2 , VARIANTEN ); IVP_STEP( X_AKT , H / 2 , VARIANTEN ); END; END; (* *) (*--------------------------------------------------------------------------*) (* *) (* IVP_HAUPT *) (* *) (*--------------------------------------------------------------------------*) (* *) BEGIN (* AWP EINLESEN *) EINLESEN; (* KOPFZEILE DRUCKEN *) WRITE_HEAD; X_AKT := 0; VRNULL := 0; VINULL := 0; MRID := ID( MRID ); (* ANFANGSWERTE FšR SPEZIELLE AUSWERTUNGEN VORBESETZEN *) MIA := MRID; MRB := MRID; MRB2 := MRID; MRB3 := MRID; MRB5 := MRID; IF AW_REAL = 0 THEN BEGIN VR_PART1 := VRY; VIY1 := 0; END ELSE BEGIN VR_PART1 := MID( VIY ); VIY1 := VIY - VR_PART1; END; VIY2 := VIY; VIY3 := VIY; VR_PART4 := VR_PART1; VIY4 := VIY1; VR_PART5 := VR_PART1; VIY5 := VIY1; VIU2 := VIY; VIU3 := VIY; VIU4 := VIY4; VIU5 := VIY5; INIT_TIMER; FAK[0] := 1; FOR I := 1 TO N DO FAK[I] := FAK[I-1] * I; FOR I := 0 TO N DO FOR J := 0 TO I DO BEGIN PROD[I-J,J] := FAK[I] / FAK[I-J]; IPROD[I-J,J] := PROD[I-J,J]; END; PRODINDEX := 0; SUBINTERVALS := 1; IVP_STEP( X_AKT , X_END - X_AKT , VARIANTEN ); CPU_TIME := GET_TIMER / 1000; WRITELN; WRITELN( 'Integrationsbereich unterteilt in ', SUBINTERVALS : 1 , ' Teilintervalle.' ); WRITELN; WRITELN( 'Rechenzeit: ', CPU_TIME : 7 : 2, ' Sekunden.' ); END; (* *) (*--------------------------------------------------------------------------*) (* *) (* HAUPTPROGRAM *) (* *) (*--------------------------------------------------------------------------*) (* *) BEGIN (* AWP EINLESEN *) EINLESEN_MAIN; IVP_HAUPT( N ); END. (*==========================================================================*) (* *) (* PASCAL-XSC - PROGRAMM PIVP9906 *) (* *) (* PROGRAMM ZUR L™SUNG DES ANFANGSWERTPROBLEMS ZUR DGL *) (* *) (* y^(n) = \sum_(i=0)^(n-1) p_i(x) y^(i) + p(x), p(x), p_i(x) Polynome *) (* *) (* POTENZREIHENANSATZ: y(x) = \sum_0^infty a_k x^k *) (* L™SUNG MIT VARIABLER TAYLORENTWICKLUNG UND RESTGLIEDABSCHŽTZUNG *) (* RECHNUNG MIT ADAPTIVER GENAUIGKEIT *) (* *) (* HIER : HAUPTPROGRAMM *) (* EINLESEN DES ANFANGSWERTPROBLEMS *) (* ZERLEGEN DES DEFINITIONSINTERVALLS IN GEEIGNETE ABSCHNITTE *) (* AUFRUF DER PROZEDUR IVP_POLY IM MODUL PIVP_MPS *) (* *) (*==========================================================================*)