ltramisc.c

linux平台下类似著名的电路板作图软件 Spice的源代码

			    *(instance->LTRAi2 + auxindex - 1),
			    *(instance->LTRAi2 + auxindex),
			    *(instance->LTRAi2 + auxindex + 1)) ;

		    eq1LTE += FABS(dashdash * h2TfirstCoeff);
		    
	    }

	    /* end LTE for convolution with i1 */

	    break;

    case LTRA_MOD_RC:

	    hilimit1 = curtime - *(ckt->CKTtimePoints + ckt->CKTtimeIndex);
	    lolimit1 = 0.0;

	    hivalue1 = LTRArcH1dashTwiceIntFunc(hilimit1,model->LTRAcByR);
	    lovalue1 = 0.0;

	    f1i = hivalue1;
	    g1i = intlinfunc(lolimit1,hilimit1,lovalue1,hivalue1,lolimit1,hilimit1);
	    h1dashTfirstCoeff = 0.5*f1i*(curtime- *(ckt->CKTtimePoints+ckt->CKTtimeIndex)) - g1i;




	    hivalue1 = LTRArcH2TwiceIntFunc(hilimit1,model->LTRArclsqr);
	    lovalue1 = 0.0;

	    f1i = hivalue1;
	    g1i = intlinfunc(lolimit1,hilimit1,lovalue1,hivalue1,lolimit1,hilimit1);
	    h1dashTfirstCoeff = 0.5*f1i*(curtime- *(ckt->CKTtimePoints+ckt->CKTtimeIndex)) - g1i;





	    hivalue1 = LTRArcH2TwiceIntFunc(hilimit1,model->LTRArclsqr);
	    lovalue1 = 0.0;

	    f1i = hivalue1;
	    g1i = intlinfunc(lolimit1,hilimit1,lovalue1,
		    hivalue1,lolimit1,hilimit1);
	    h1dashTfirstCoeff = 0.5*f1i*(curtime-
		*(ckt->CKTtimePoints+ckt->CKTtimeIndex)) - g1i;



	    /* LTEs for convolution with v1 */
	    /* get divided differences for v1 (2nd derivative estimates) */

	    /* 
	     * no need to subtract operating point values because
	     * taking differences anyway
	     */

	    dashdash = SECONDDERIV(ckt->CKTtimeIndex+1,
			    *(instance->LTRAv1+ckt->CKTtimeIndex-1),
			    *(instance->LTRAv1+ckt->CKTtimeIndex),
			    *(ckt->CKTrhsOld + instance->LTRAposNode1) -
			    *(ckt->CKTrhsOld + instance->LTRAnegNode1));
	    eq1LTE += FABS(dashdash * h1dashTfirstCoeff);
	    eq2LTE += FABS(dashdash * h3dashTfirstCoeff);
		    
	    /* end LTEs for convolution with v1 */

	    /* LTEs for convolution with v2 */
	    /* get divided differences for v2 (2nd derivative estimates) */

	    dashdash = SECONDDERIV(ckt->CKTtimeIndex+1,
			    *(instance->LTRAv2+ckt->CKTtimeIndex-1),
			    *(instance->LTRAv2+ckt->CKTtimeIndex),
			    *(ckt->CKTrhsOld + instance->LTRAposNode2) -
			    *(ckt->CKTrhsOld + instance->LTRAnegNode2));

	    eq2LTE += FABS(dashdash * h1dashTfirstCoeff);
	    eq1LTE += FABS(dashdash * h3dashTfirstCoeff);
		    
	    /* end LTEs for convolution with v2 */

	    /* LTE for convolution with i1 */
	    /* get divided differences for i1 (2nd derivative estimates) */

	    dashdash = SECONDDERIV(ckt->CKTtimeIndex+1,
				    *(instance->LTRAi1 + ckt->CKTtimeIndex - 1),
				    *(instance->LTRAi1 + ckt->CKTtimeIndex),
				    *(ckt->CKTrhsOld + instance->LTRAbrEq1));

	    eq2LTE += FABS(dashdash * h2TfirstCoeff);
		    

	    /* end LTE for convolution with i1 */

	    /* LTE for convolution with i2 */
	    /* get divided differences for i2 (2nd derivative estimates) */

	    dashdash = SECONDDERIV(ckt->CKTtimeIndex+1,
					    *(instance->LTRAi2 + ckt->CKTtimeIndex - 1),
					    *(instance->LTRAi2 + ckt->CKTtimeIndex),
					    *(ckt->CKTrhsOld + instance->LTRAbrEq2));


	    eq1LTE += FABS(dashdash * h2TfirstCoeff);
		    

	    /* end LTE for convolution with i1 */

	    break;

    default: return(1/*error*/);
    }


#ifdef LTRADEBUG
	fprintf(stdout,"%s: LTE/input for Eq1 at time %g is: %g\n",
		instance->LTRAname, curtime, eq1LTE/instance->LTRAinput1);
	fprintf(stdout,"%s: LTE/input for Eq2 at time %g is: %g\n",
		instance->LTRAname, curtime, eq2LTE/instance->LTRAinput1);
	fprintf(stdout,"\n");
#endif

    return(FABS(eq1LTE) + FABS(eq2LTE));
}

/*********************************************************************/
/****************** old stuff, retained for historical interest ******/
/*********************************************************************/



/*
 * LTRAcoeffSetup sets up the coefficient list for the convolution,
 * returns the coefficient at (current_timepoint-T)
 */

/*
double 
LTRAcoeffSetup(coefflist,listsize,T,firstvalue,valuelist,curtime,timelist,timeindex,auxindexptr)
double *coefflist, *timelist, *valuelist;
int listsize, timeindex;
double T, firstvalue, curtime;
int *auxindexptr;

{
unsigned exact;
double returnval, delta1, delta2;
double dummy1, dummy2;
double lolimit1,lolimit2,hilimit1,hilimit2;
double lovalue1,lovalue2,hivalue1,hivalue2;
int i,auxindex;
*/

/* coefflist should already have been allocated to the necessary size */

/*
#ifdef LTRAdebug
if (listsize <= timeindex) {
	printf("LTRAcoeffSetup: not enough space in coefflist\n");
}
#endif

*/

/*
 * we assume a piecewise linear function, and we calculate the
 * coefficients using this assumption in the integration of the
 * function
 */

/*
if (T == 0.0) {
	auxindex = timeindex;
} else {

	if (curtime - T <= 0.0) {
		for (i =0; i<= timeindex; i++) {
			*(coefflist + i) = 0.0;
		}
		*auxindexptr = 0;
		return(0.0);
	} else {
		exact = 0;
		for (i = timeindex; i>= 0; i--) {
			if (curtime - *(timelist + i) ==  T) {
				exact =1;
				break;
			} 
			if (curtime - *(timelist + i) > T) break;
		}

#ifdef LTRADEBUG
		if ((i < 0) || ((i==0) && (exact==1)))
			printf("LTRAcoeffSetup: i <= 0: some mistake!\n");
#endif

		if (exact == 1) {
			auxindex = i-1;
		} else {
			auxindex = i;
		}
	}
}
*/
/* the first coefficient */

/*
delta1 = curtime -T - *(timelist + auxindex);
lolimit1 = T;
hilimit1 = T + delta1;
lovalue1 = firstvalue;
hivalue1 = *(valuelist + auxindex);
dummy1 = twiceintlinfunc(lolimit1,hilimit1,lolimit1,lovalue1,
	hivalue1,lolimit1,hilimit1)/delta1;
returnval = dummy1;

*/


/* the coefficients for the rest of the timepoints */

/*
for (i=auxindex; i>0; i--) {

	delta2 = delta1;*/ /* previous delta1 */
	/*lolimit2 = lolimit1;*/ /* previous lolimit1 */
	/*hilimit2 = hilimit1;*/ /*previous hilimit1 */
	/*lovalue2 = lovalue1;*/ /* previous lovalue1 */
	/*hivalue2 = hivalue1;*/ /*previous hivalue1 */
	/*dummy2 = dummy1;*/ /* previous dummy1 */
	/*delta1 = *(timelist + i) - *(timelist + i - 1);
	lolimit1 = hilimit2;
	hilimit1 = curtime - *(timelist + i - 1);
	lovalue1 = hivalue2;
	hivalue1 = *(valuelist + i - 1);
	dummy1 = twiceintlinfunc(lolimit1,hilimit1,lolimit1,
		lovalue1,hivalue1,lolimit1,hilimit1)/delta1;

	*(coefflist + i) = dummy1 - dummy2 + intlinfunc(lolimit2,hilimit2,
		lovalue2,hivalue2,lolimit2,hilimit2);
}
*auxindexptr = auxindex;
return(returnval);
}
*/

/*
 * LTRAtCoeffSetup sets up the coefficient list for the LTE calculation,
 * returns the coefficient at (current_timepoint-T)
 */

/*
double 
LTRAtCoeffSetup(coefflist,listsize,T,valuelist,
	firstothervalue,othervaluelist,curtime,timelist,timeindex,
	auxindexptr, ltecontype)
double *coefflist, *timelist, *valuelist, *othervaluelist;
int listsize, timeindex;
double T, firstothervalue, curtime;
int *auxindexptr, ltecontype;

{
unsigned exact;
double returnval, delta;
double dummy;
double f1i, f2i, g1i, g2i;
double lolimit1, hilimit1;
double lovalue1, hivalue1;
double lolimit2, hilimit2;
double lovalue2, hivalue2;
double firstint1 = 0.0, firstint2 = 0.0;
double secondint1 = 0.0, secondint2 = 0.0;
int i,auxindex;

*/
/* coefflist should already have been allocated to the necessary size */

/*
#ifdef LTRAdebug
if (listsize <= timeindex) {
	printf("LTRAtCoeffSetup: not enough space in coefflist\n");
}
#endif

*/

/*
 * we assume a piecewise linear function, and we calculate the
 * coefficients using this assumption in the integration of the
 * function
 */

/*
if (T == 0.0) {
	auxindex = timeindex;
} else {

	if (curtime - T <= 0.0) {
		for (i =0; i<= timeindex; i++) {
			*(coefflist + i) = 0.0;
		}
		*auxindexptr = 0;
		return(0.0);
	} else {
		exact = 0;
		for (i = timeindex; i>= 0; i--) {
			if (curtime - *(timelist + i) ==  T) {
				exact =1;
				break;
			} 
			if (curtime - *(timelist + i) > T) break;
		}

#ifdef LTRADEBUG
		if ((i < 0) || ((i==0) && (exact==1)))
			printf("LTRAcoeffSetup: i <= 0: some mistake!\n");
#endif

		if (exact == 1) {
			auxindex = i-1;
		} else {
			auxindex = i;
		}
	}
}
*/
/* the first coefficient */

/* i = n in the write-up */

/*
hilimit1 = curtime - *(timelist + auxindex);
hivalue1 = *(valuelist + auxindex);
lolimit1 = *(timelist + timeindex) - *(timelist + auxindex);
lolimit1 = MAX(T,lolimit1);
lovalue1 = firstothervalue;
f1i = twiceintlinfunc(lolimit1,hilimit1,lolimit1,lovalue1,hivalue1,lolimit1,
	hilimit1);
g1i = thriceintlinfunc(lolimit1,hilimit1,lolimit1,lolimit1,lovalue1,
	hivalue1,lolimit1,hilimit1);
returnval = 0.5*f1i*(curtime-T- *(timelist+auxindex)) - g1i;

*/
/* the coefficients for the rest of the timepoints */

/*
if (ltecontype != LTRA_MOD_HALFCONTROL) {
for (i=auxindex; i>0; i--) {

	lolimit2 = lolimit1;*/ /* previous lolimit1 */
	/*hilimit2 = hilimit1;*/ /*previous hilimit1 */
	/*lovalue2 = lovalue1;*/ /* previous lovalue1 */
	/*hivalue2 = hivalue1;*/ /*previous hivalue1 */
	/*f2i = f1i;*/ /* previous f1i */
	/*g2i = g1i;*/ /* previous g1i */
	/*firstint2 = firstint1;*/ /* previous firstint1 */
	/*secondint2 = secondint1;*/ /* previous secondint1 */
	/*lolimit1 = *(timelist + timeindex) - *(timelist + i - 1);
	hilimit1 = curtime - *(timelist + i - 1);
	lovalue1 = *(othervaluelist + i - 1);
	hivalue1 = *(valuelist + i - 1);
	firstint1 += intlinfunc(lolimit2, lolimit1, lovalue2, lovalue1,
		lolimit2, lolimit1);
	secondint1 += (lolimit1-lolimit2)*firstint2 + twiceintlinfunc(
		lolimit2,lolimit1,lolimit2,lovalue2,lovalue1,lolimit2,lolimit1);
	f1i = twiceintlinfunc(lolimit1,hilimit1,lolimit1,lovalue1,hivalue1,
		lolimit1,hilimit1) + firstint1*(hilimit1-lolimit1);
	g1i = thriceintlinfunc(lolimit1,hilimit1,lolimit1,lolimit1,
		lovalue1,hivalue1,lolimit1,hilimit1) +
		(hilimit1-lolimit1)*(hilimit1-lolimit1)*0.5*firstint1 +
		(hilimit1-lolimit1)*secondint1;

	*(coefflist + i) = g2i - g1i + 0.5*(f1i + f2i)*(*(timelist+i) -
		*(timelist+i-1));
}
}
*auxindexptr = auxindex;
return(returnval);
}
*/


/* formulae taken from the Handbook of Mathematical Functions by
 * Milton Abramowitz and Irene A. Stegan, page 378, formulae 9.8.1 -
 * 9.8.4
 */

/*

double
bessi0(x)
double x;
{
	double t, tsq, oneovert, result, dummy;
	int i;
	static double coeffs1[7], coeffs2[9];

	coeffs1[0] = 1.0;
	coeffs1[1] = 3.5156229;
	coeffs1[2] = 3.0899424;
	coeffs1[3] = 1.2067492;
	coeffs1[4] = 0.2659732;
	coeffs1[5] = 0.0360768;
	coeffs1[6] = 0.0045813;

	coeffs2[0] = 0.39894228;
	coeffs2[1] = 0.01328592;
	coeffs2[2] = 0.00225319;
	coeffs2[3] = -0.00157565;
	coeffs2[4] = 0.00916281;
	coeffs2[5] = -0.02057706;
	coeffs2[6] = 0.02635537;
	coeffs2[7] = -0.01647633;
	coeffs2[8] = 0.00392377;

	t = x/3.75;
	dummy = 1.0;

	if (fabs(t) <= 1) {
		tsq = t*t;

		result = 1.0;
		for (i=1;i<=6;i++) {
			dummy *= tsq;
			;
			result += dummy * coeffs1[i];
		}
	}
	else
	{
		oneovert = 1/fabs(t);

		result = coeffs2[0];
		for (i=1;i<=8;i++) {
			dummy *= oneovert;
			result += coeffs2[2] * dummy;
		}
		result *= exp(x) * sqrt(1/fabs(x));
	}
	return(result);
}

double
bessi1(x)
double x;
{
	double t, tsq, oneovert, result, dummy;
	int i;
	static double coeffs1[7], coeffs2[9];

	coeffs1[0] = 0.5;
	coeffs1[1] = 0.87890594;
	coeffs1[2] = 0.51498869;
	coeffs1[3] = 0.15084934;
	coeffs1[4] = 0.02658733;
	coeffs1[5] = 0.00301532;
	coeffs1[6] = 0.00032411;

	coeffs2[0] = 0.39894228;
	coeffs2[1] = -0.03988024;
	coeffs2[2] = -0.00362018;
	coeffs2[3] = 0.00163801;
	coeffs2[4] = -0.01031555;
	coeffs2[5] = 0.02282967;
	coeffs2[6] = -0.02895312;
	coeffs2[7] = 0.01787654;
	coeffs2[8] = -0.00420059;

	t = x/3.75;
	dummy = 1.0;

	if (fabs(t) <= 1) {
		tsq = t*t;

		result = 0.5;
		for (i=1;i<=6;i++) {
			dummy *= tsq;
			;
			result += dummy * coeffs1[i];
		}
		result *= x;
	}
	else
	{
		oneovert = 1/fabs(t);

		result = coeffs2[0];
		for (i=1;i<=8;i++) {
			dummy *= oneovert;
			result += coeffs2[2] * dummy;
		}
		result *= exp(x) * sqrt(1/fabs(x));
		if (x < 0) result = -result;
	}
	return(result);
}
*/

/*
 * LTRAdivDiffs returns divided differences after 2 iterations,
 * an approximation to the second derivatives. The algorithm is
 * picked up directly from Tom Quarles' CKTterr.c; no attempt has
 * been made to figure out why it does what it does.
 */

/*
double
LTRAdivDiffs(difflist, valuelist, firstvalue, curtime, timelist, 
	timeindex)
double *difflist, *valuelist, firstvalue,  *timelist, curtime;
int timeindex;

{
double *dtime, *diffs, returnval;
int i,j;