WilsonWebs - An automated tool for n-Loop Wilson web generation

git clone git://git.bcharge.de/WilsonWebs.git

About | Log | Files | Refs | License

webs.c (10420B)


/* See LICENSE file for copyright and license details. */
#include <stdio.h>
#include <string.h>

#include "util.h"

typedef struct {
	unsigned int legId;
	unsigned int id;
} Vertex;

typedef struct {
	char *name;
	int init;
	unsigned int id;
	int isMassive;
	int isGluon;
	Vertex *Vertices;
} Leg;


typedef struct {
	Vertex emitter;
	Vertex absorber;
} Link;

static int nGraphs = 0;

void
generate_legs(char **init, size_t n_init, char **fin, size_t n_fin, size_t nLoops, Leg *Legs)
{
	for (int i=0; i<n_init; ++i) {
		int isGluon = 0;
		if (!strcmp("g", init[i]))
			isGluon = 1;
		Vertex *Vertices = (Vertex*)ecalloc(nLoops+2, sizeof(Vertex));
		Vertices[0] = (Vertex){i+1,0};
		Legs[i+1] = (Leg){init[i], 1, i, 0, isGluon, Vertices};
	}

	for (int i=0; i<n_fin; ++i) {
		int isGluon = 0;
		if (!strcmp("g", fin[i]))
			isGluon = 1;
		Vertex *Vertices = (Vertex *)ecalloc(nLoops+2, sizeof(Vertex));
		Vertices[0] = (Vertex){n_init+i+1,0};
		Legs[n_init+i+1] = (Leg){fin[i], 0, n_init+i, 0, isGluon,Vertices};
	}
}

void
generate_vertices(Leg *Legs, size_t nLegs, int nLoops, Vertex *Vertices)
{
	for (int i=0; i<nLegs; ++i) {
		for (int j=0; j<nLoops; ++j) {
			/* TODO: Don't forget to add vertices for self-energies later */
			Vertices[j+nLoops*i] = (Vertex){i+1, 1+j};
			Legs[i+1].Vertices[j+1] = Vertices[j+nLoops*i];
		}
			Vertex V = (Vertex){i+1,999};
			Legs[i+1].Vertices[nLoops+1] = V;
	}

	/*
	for (int i = 0; i < nLegs*nLoops; ++i) {
		printf("test %d\n", Vertices[i].legId);
	}
	*/
}

void
generate_loops(Vertex *Vertices, Leg *Legs, Link **Graphs, size_t nVx, size_t nLegs, int nLoops, int nGenLoops)
{
	if (nGenLoops == nLoops) {
		printf("End recursion.\n");
		free(Vertices);
		nGraphs++;
		return;
	}
	for (int iVx = 0; iVx < nVx; ++iVx) {
		Vertex V = Vertices[iVx];
		size_t nTargets = nLegs-V.legId;
		Vertex *targets = (Vertex*)ecalloc(nTargets, sizeof(Vertex));
		int iTar = -1;
		for (int iLeg = V.legId+1; iLeg <= nLegs; ++iLeg) {
			/*printf("Leg %d\n", iLeg);*/
			targets[++iTar] = Legs[iLeg].Vertices[V.id];
		}

		for (int iTar = 0; iTar < nTargets; ++iTar) {
			/*printf("cycle\n");*/
			Vertex T = targets[iTar];
			printf("Loop: %d\tV%d%d - V%d%d\n", nGenLoops, V.legId, V.id, T.legId, T.id);
			Graphs[nGraphs][nGenLoops] = (Link){V,T};
			/*
			int i = 0;
			int iGraph = nGraphs;
			printf("Graph# %d\tV%d%d - V%d%d\n", iGraph, Graphs[iGraph][i].emitter.legId, Graphs[iGraph][i].emitter.id, Graphs[iGraph][i].absorber.legId, Graphs[iGraph][i].absorber.id);
			*/
			//targets += 1;
			nGenLoops++;
			printf("nGenLoops %d\n", nGenLoops);
			//generate_loops(targets, Legs, Graphs, nTargets-1, nLegs, nLoops, nGenLoops--);
			Vertex *remainingVtcs = (Vertex *)ecalloc(nVx-1, sizeof(Vertex));
			memmove(&remainingVtcs[0], &Vertices[1], (nVx - 1) * sizeof(Vertex));
			generate_loops(remainingVtcs, Legs, Graphs, nTargets-1, nLegs, nLoops, nGenLoops--);
		}
	}
	//generate_loops(Vertices, deletedVertices, Graphs, nVx, nLoops, nGenLoops);
	////memmove(&AllPossVertices[0], &Vertices[1], (n_vx - 1) * sizeof(Vertex));

}


int
main(int argc, char *argv[])
{
	if (argc == 2 && !strcmp("-v", argv[1]))
		die("WilsonWebs-"VERSION);
	else if (argc != 1)
		die("usage: webs [-v]");

	char *initial[] = {"q","qbar"};
	//char *initial[] = {"q","qbar"};
	//char *final[]   = {"q","qbar"};
	char *final[]   = {"q","qbar", "q", "qbar"};
	//char *final[]   = {"q","qbar", "q", "qbar", "q", "qbar","q", "qbar","q", "qbar","q","qbar", "q", "qbar", "q", "qbar","q", "qbar"};

	FILE *fcol = fopen("colorFactors.h","w");
	FILE *fkin = fopen("kinFactors.h","w");

	if (fcol == NULL || fkin == NULL) {
		die("Error while creation files for writing!");
	}

	int nLoops = 1;

	size_t nLegs = NELEMS(initial) + NELEMS(final);
	Leg Legs[nLegs+1];
	generate_legs(initial,NELEMS(initial), final, NELEMS(final), nLoops, Legs);

	size_t nVx = nLegs*nLoops;
	Vertex *Vertices = (Vertex *)ecalloc(nVx, sizeof(Vertex));
	generate_vertices(Legs, nLegs, nLoops, Vertices);

	for (int i = 0; i < nVx; ++i) {
		printf("V%d%d\n",Vertices[i].legId,Vertices[i].id);
	}

	/*int *deletedVertices = (int *)ecalloc(nVx, sizeof(int));*/
	Link **Graphs = (Link **)ecalloc(300*nLoops, sizeof(Link*));
	for (int i = 0; i < 300; ++i) {
		Graphs[i] = (Link*)ecalloc(nLoops, sizeof(Link));
	}

	Link *Skeleton = (Link *)ecalloc(300*nLoops, sizeof(Link));
	int iEdge = -1;
	for (int iLeg = 1; iLeg <= nLegs; ++iLeg) {
		Vertex* LegVertices = Legs[iLeg].Vertices;
		Skeleton[++iEdge] = (Link){LegVertices[0],LegVertices[1]};
		for (int iVx = 1; iVx <= nLoops; ++iVx) {
			Skeleton[++iEdge] = (Link){LegVertices[iVx],LegVertices[iVx+1]};
		}
		/*printf("Last vx: V%d%d\n", LegVertices[nLoops].legId, LegVertices[nLoops].id);*/
		//Skeleton[++iEdge] = (Link){LegVertices[nLoops], (Vertex){999, 999}};
	}

	for (int i = 0; i <= iEdge; ++i) {
		printf("Parton lines: Emit V%d%d Absorb V%d%d\n",
				Skeleton[i].emitter.legId,Skeleton[i].emitter.id, Skeleton[i].absorber.legId,Skeleton[i].absorber.id);
	}

	printf("Edge count: %d\n", iEdge);

	generate_loops(Vertices, Legs, Graphs, nVx, nLegs, nLoops, 0);
	for (int iGraph=0; iGraph<nGraphs; ++iGraph) {
		for (int i = 0; i < nLoops; ++i) {
			printf("Graph# %d\tV%dx%d - V%dx%d\n", iGraph, Graphs[iGraph][i].emitter.legId, Graphs[iGraph][i].emitter.id, Graphs[iGraph][i].absorber.legId, Graphs[iGraph][i].absorber.id);
			int emitterLegId = Graphs[iGraph][i].emitter.legId;
			int emitterId = Graphs[iGraph][i].emitter.id;
			Vertex *emitterVtcs = Legs[emitterLegId].Vertices;
			int emitterNxtId = emitterVtcs[emitterId+1].id;

			int absorberLegId = Graphs[iGraph][i].absorber.legId;
			int absorberId = Graphs[iGraph][i].absorber.id;
			Vertex *absorberVtcs = Legs[absorberLegId].Vertices;
			int absorberNxtId = absorberVtcs[absorberId+1].id;

			char *color = (char*)malloc(10000 * sizeof(char));
			char *kin = (char*)malloc(10000 * sizeof(char));
			char *Delta1 = (char*)malloc(sizeof(char));
			char *Delta2 = (char*)malloc(sizeof(char));
			char *delta1 = (char*)malloc(sizeof(char));
			char *delta2 = (char*)malloc(sizeof(char));

			if(!strcmp(Legs[emitterLegId].name,"q")){
				Delta1 = "+";
			}
			else if(!strcmp(Legs[emitterLegId].name,"qbar")){
				Delta1 = "-";
			}
			if(!strcmp(Legs[absorberLegId].name,"q")){
				Delta2 = "+";
			}
			else if(!strcmp(Legs[absorberLegId].name,"qbar")){
				Delta2 = "-";
			}

			if(Legs[emitterLegId].init && Legs[absorberLegId].init)
			{
				delta1 = "-";
				delta2 = "+";
				if(!strcmp(Legs[emitterLegId].name,"g"))
					Delta1 = "-";
				if(!strcmp(Legs[absorberLegId].name,"g"))
					Delta2 = "-";
			}
			else if(Legs[emitterLegId].init && !Legs[absorberLegId].init)
			{
				delta1 = "-";
				delta2 = "-";
				if(!strcmp(Legs[emitterLegId].name,"g"))
					Delta1 = "-";
				if(!strcmp(Legs[absorberLegId].name,"g"))
					Delta2 = "-";
			}
			else if(!Legs[emitterLegId].init && !Legs[absorberLegId].init)
			{
				delta1 = "+";
				delta2 = "-";
				if(!strcmp(Legs[emitterLegId].name,"g"))
					Delta1 = "-";
				if(!strcmp(Legs[absorberLegId].name,"g"))
					Delta2 = "-";
			}

			sprintf(color, "L C%dx%d = ", emitterLegId, absorberLegId);
			sprintf(kin,   "L A%dx%d = ", emitterLegId, absorberLegId);
			sprintf(kin+strlen(kin), "(%sv%d(n1))/(%sv%d.l)*(%sv%d(m1))/(%sv%d.l)*Ngluon(n1,m1)/l.l;\n",Delta1,emitterLegId,delta1,emitterLegId,Delta2,absorberLegId,delta2,absorberLegId);

			if (Legs[emitterLegId].init && !strcmp(Legs[emitterLegId].name,"q")){
				sprintf(color+strlen(color), "SUNT(Glu99,Col%d%d,Col%d)*",emitterLegId,emitterNxtId,emitterLegId);
			}
			else if (!Legs[emitterLegId].init && !strcmp(Legs[emitterLegId].name,"q")){
				sprintf(color+strlen(color), "SUNT(Glu99,Col%d,Col%d%d)*",emitterLegId,emitterLegId,emitterNxtId);
			}
			else if (Legs[emitterLegId].init && !strcmp(Legs[emitterLegId].name,"qbar")){
				sprintf(color+strlen(color), "SUNT(Glu99,Col%d,Col%d%d)*",emitterLegId,emitterLegId,emitterNxtId);
			}
			else if (!Legs[emitterLegId].init && !strcmp(Legs[emitterLegId].name,"qbar")){
				sprintf(color+strlen(color), "SUNT(Glu99,Col%d%d,Col%d)*",emitterLegId,emitterNxtId,emitterLegId);
			}

			if (Legs[absorberLegId].init && !strcmp(Legs[absorberLegId].name,"q")){
				sprintf(color+strlen(color), "SUNT(Glu99,Col%d%d,Col%d)*",absorberLegId,absorberNxtId,absorberLegId);
			}
			else if (!Legs[absorberLegId].init && !strcmp(Legs[absorberLegId].name,"q")){
				sprintf(color+strlen(color), "SUNT(Glu99,Col%d,Col%d%d)*",absorberLegId,absorberLegId,absorberNxtId);
			}
			else if (Legs[absorberLegId].init && !strcmp(Legs[absorberLegId].name,"qbar")){
				sprintf(color+strlen(color), "SUNT(Glu99,Col%d,Col%d%d)*",absorberLegId,absorberLegId,absorberNxtId);
			}
			else if (!Legs[absorberLegId].init && !strcmp(Legs[absorberLegId].name,"qbar")){
				sprintf(color+strlen(color), "SUNT(Glu99,Col%d%d,Col%d)*",absorberLegId,absorberNxtId,absorberLegId);
			}

			if (Legs[emitterLegId].init && !strcmp(Legs[emitterLegId].name,"g")){
				sprintf(color+strlen(color), "(-i_)*SUNF(Glu99,Glu%d,Glu%d%d)*",emitterLegId,emitterLegId,emitterNxtId);
			}
			else if (!Legs[emitterLegId].init && !strcmp(Legs[emitterLegId].name,"g")){
				sprintf(color+strlen(color), "(-i_)*SUNF(Glu99,Glu%d%d,Glu%d)*",emitterLegId,emitterNxtId,emitterLegId);
			}

			if (Legs[absorberLegId].init && !strcmp(Legs[absorberLegId].name,"g")){
				sprintf(color+strlen(color), "(-i_)*SUNF(Glu99,Glu%d,Glu%d%d)*",absorberLegId,absorberLegId,absorberNxtId);
			}
			else if (!Legs[absorberLegId].init && !strcmp(Legs[absorberLegId].name,"g")){
				sprintf(color+strlen(color), "(-i_)*SUNF(Glu99,Glu%d%d,Glu%d)*",absorberLegId,absorberNxtId,absorberLegId);
			}

			for (int iLeg = 1; iLeg <= nLegs; ++iLeg) {
				if (iLeg != emitterLegId && iLeg !=absorberLegId) {
					sprintf(color+strlen(color), "d_(Col%d,Col%d%d)*",Legs[iLeg].Vertices[0].legId,Legs[iLeg].Vertices[0].legId,Legs[iLeg].Vertices[nLoops+1].id );
				}
			}
			sprintf(color+strlen(color), "1;\n");
			fprintf(fcol,"%s", color);
			fprintf(fkin,"%s", kin);
		}
	}
	fclose(fcol);
	fclose(fkin);
	return 0;

	/*
	for (int iVx=0; iVx<300; ++iVx) {
		printf("Graph ID: %d\n",iVx);
		for (int i=0; i<nLoops; ++i) {
			printf("V%d%d-V%d%d\n",Graphs[iVx*nLoops+i].emitter.legId, Graphs[iVx*nLoops+i].emitter.id, Graphs[iVx*nLoops+i].absorber.legId, Graphs[iVx*nLoops+i].absorber.id);
		}
	}
	*/

	//draw(Legs, Connections);


	return 0;
}