[..]

Quote:

>>> On my old P54c 166 MHz machine your benchmark (the gForth

>>> variant) executes in 1.167 seconds. This should be roughly

>>> equivalent to 587 ms on a 330 MHz PIII, [...]

[..]

Here is slbench with routines from the FSL and a plotter-type

visualization. I'm using the variable step Runge-Kutta ODE solver

with an eps of 1e-3. Instead of 1.167 seconds it now takes 15 ms

(on my P54c) to get the (nice-looking) output.

-marcel

-- -----------

ANEW -slbench3

NEEDS -rk4

NEEDS -sketcher

DOC

(*

Benchmark version of sl for iForth

Solve the semiconductor laser rate equations, given a current pulse

profile. Output the time, intensity, phase, and current density.

Krishna Myneni, 1-26-2000

*)

ENDDOC

\ ================

\ Laser parameters

\ ================

4.5e-12 FVALUE t_p \ photon lifetime (sec)

700e-12 FVALUE t_s \ carrier lifetime (sec)

2.6e-6 FVALUE G_N \ differential gain (cm^3/s)

1.5e18 FVALUE N_th \ threshold carrier density (cm^-3)

20e-3 FVALUE I_th \ threshold current through laser (mA)

5e FVALUE alpha \ linewidth enhancement factor

\ ========================

\ Dimensionless parameters

\ ========================

0e FVALUE T_ratio \ T_ratio = t_s/t_p

0e FVALUE PumpFactor \ P = PumpFactor*(I/I_th - 1)

0e FVALUE P

\ compute the normalized parameters

: init_params ( -- )

t_s t_p F/ TO T_ratio

t_p G_N F* N_th F* F2/ TO PumpFactor ;

\ =============================

\ The injection current profile

\ =============================

1e-9 FVALUE fwhm \ full-width at half-max for

\ current pulse in ns (1 ns)

20e-3 FVALUE pulse_amp \ current pulse amplitude

\ above d.c. level (20 mA)

I_th 10e-3 F+ FVALUE dc_current \ d.c. current level (10 mA

\ above threshold)

50e-9 FVALUE peak_offset \ offset for current peak (50 ns)

\ compute current at real time t

:INLINE GaussianPulse ( F: t -- c )

peak_offset F-

fwhm F/

FSQR -2.77066e F* FEXP

pulse_amp F* dc_current F+ ;

\ =============================

\ Pump rate

\ =============================

\ P = PumpFactor*(I/I_th - 1);

:INLINE pump ( F: c -- )

I_th F/ F1-

PumpFactor F* TO P ;

\ ======================================

\ Rates of dimensionless state variables

\ ======================================

\ Data in a state vector has the following order: Re{Y}, Im{Y}, Z

0 VALUE 's

0 VALUE 't

:INLINE t[]! 't []DFLOAT DF! ;

:INLINE intensity^2 ( -- ) ( F: -- I ) \ compute intensity

0 s[] FSQR 1 s[] FSQR F+ ;

:INLINE intensity ( -- ) ( F: -- I ) \ compute intensity

intensity^2 FSQRT ;

:INLINE phase ( -- ) ( F: -- phase ) \ compute phase in radians

1 s[] 0 s[] FATAN2 ;

\ compute dY/ds for state vector 'a'

\ Re{dY/ds} = Z(Re{Y} - alpha*Im{Y})

\ Im{dY/ds} = Z(alpha*Re{Y} + Im{Y})

: dY/ds ( F: -- )

0 s[] 1 s[] alpha F* F- 2 s[] F* 0 t[]!

0 s[] alpha F* 1 s[] F+ 2 s[] F* 1 t[]! ;

\ compute rate of change of Z

\ dZ/ds = (P - Z - (1 + 2Z)|Y|^2)/T

: dZ/ds ( F: -- )

intensity^2 2 s[] F2* F1+ F*

P 2 s[] F- FSWAP F- T_ratio F/ 2 t[]! ;

\ derivative of the state vector

:NONAME ( u{ dudt{ -- ) ( F: t -- )

TO 't TO 's GaussianPulse pump dY/ds dZ/ds ; =: 'dsdt

\ ==========

\ ODE Solver

\ ==========

1e FVALUE ds

25000e t_p F* 1e-9 F/ FVALUE endtime

0 VALUE ix

3 DOUBLE ARRAY V{

#512 DOUBLE ARRAY t{

#512 DOUBLE ARRAY i{ -- intensity

#512 DOUBLE ARRAY p{ -- phase

#512 DOUBLE ARRAY c{ -- current density

: main ( -- )

CLEAR ix

init_params

2e FSQRT V{ 0 } DF!

0e V{ 1 } DF!

0e V{ 2 } DF!

'dsdt 3 V{ ds 1e-3 )rk4qc_init

ds 0e

BEGIN

rk4qc_step 0= IF CR ." Time: " F.N1

." s, step: " F.N1

." s, ix = " ix DEC.

TRUE ABORT" Time step too small"

ENDIF

ix #512 < IF FDUP t{ ix } DF!

intensity i{ ix } DF!

phase p{ ix } DF!

2 s[] c{ ix } DF!

ENDIF

1 +TO ix FDUP endtime F>

UNTIL

F2DROP rk4qc_done

CR ." needed " ix DEC. ." steps. " ;

CR TIMER-RESET main .ELAPSED

\ EOF