#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <cctype>
#include <cmath>
#include <algorithm>
#define rep(i, a, b) for (int i = (a), _ = (b); i <= _; ++ i)
#define per(i, a, b) for (int i = (a), _ = (b); i >= _; -- i)
#define For(i, a, b) for (int i = (a), _ = (b); i < _; ++ i)
#define fore(e, u, v) for (int p = e(u), v = e[p].y; p; v = e[p = e[p].nxt].y)
#define ri rd<int>
using namespace std;
const int maxN = 1e5 + 7;
const int mod = 1e9 + 7;

inline int pls(int x, int y) {return (x + y) % mod;}
inline int mul(int x, int y) {return 1LL * x * y % mod;}
inline int mns(int x, int y) {return pls(x, mod - y);}
inline void Mul(int &x, int y) {x = mul(x, y);}
inline void Add(int &x, int y) {x = pls(x, y);}

template<class T> inline T rd() {
	bool f = 1; char c = getchar(); for (; !isdigit(c); c = getchar()) if (c == '-') f = 0;
	T x = 0; for (; isdigit(c); c = getchar()) x = x * 10 + c - 48; return f ? x : -x;
}

inline void NIE() {puts("0"); exit(0);}

int n;
int c[maxN], to[maxN];
int ans = 1;

namespace G {

	int vis[maxN];
	int lop[maxN], tl;
	int stack[maxN], Top;
	int len[maxN], s[maxN];

	struct vec {
		static const int maxE = maxN;
		int g[maxN], te;
		struct edge {int y, nxt, fl;} e[maxE << 1];
		vec() {te = 1;}
		inline void push(int x, int y) {e[++te] = (edge){y, g[x], 0}; g[x] = te;}
		inline void link(int x, int y) {push(x, y), push(y, x);}
		inline int& operator () (int x) {return g[x];}
		inline edge& operator [] (int x) {return e[x];}
	}e;

	bool dfs(int x) {
		stack[++Top] = x;
		vis[x] = 1;
		fore (e, x, y) if (e[p].fl == 0) {
			e[p].fl = e[p ^ 1].fl = 1;
			if (!vis[y]) {
				if (dfs(y)) return 1;
			}
			else {
				do {
					vis[stack[Top]] = 2;
					lop[++tl] = stack[Top];
				} while (stack[Top--] != y);
				return 1;
			}
		}
		stack[Top--] = 0;
		return 0;
	}

	int depth(int x, int fa) {
		int d = 1;
		if (fa) vis[x] = 1;
		fore (e, x, y) if (y != fa && vis[y] != 2) {
			if (d != 1) NIE();
			d = 1 + depth(y, x);
		}
		return d;
	}

	void gao() {
		rep (x, 1, n) if (!vis[x]) {
			tl = Top = 0;
			dfs(x);

			if (tl > 1 && to[lop[2]] == lop[1]) reverse(lop+1, lop+tl+1);
			rep (i, 1, tl) len[i] = len[i + tl] = depth(lop[i], 0) - 1;
			rep (i, 1, tl) if (len[i] > tl) NIE();

			int cnt = 0;
			rep (i, 1, tl) if (len[i] != 0) ++ cnt;
			if (cnt == 0) ++c[tl]; 
			else {
				for (int i = tl << 1; i > tl; -- i) if (len[i]) {
					int l = i - len[i]; 
					for (-- i; i > l; -- i) if (len[i]) NIE();
					if (len[i]) Mul(ans, 1), ++ i;
					else Mul(ans, 2);
				}
			}
		}
	}
}

namespace L {
	int f[maxN], fac[maxN], ifac[maxN], pw2[maxN], iv[maxN];

	void init() {
		iv[1] = 1; rep (i, 2, n) iv[i] = mul(iv[mod % i], mod - mod / i);
		fac[0] = 1; rep (i, 1, n) fac[i] = mul(fac[i - 1], i);
		ifac[0] = 1; rep (i, 1, n) ifac[i] = mul(ifac[i - 1], iv[i]);
		pw2[0] = 1; rep (i, 1, n) pw2[i] = mul(pw2[i - 1], 2);
	}

	void solve() {
		init();

		int res = 1;

		f[0] = f[1] = 1;
		rep (i, 2, n) f[i] = pls(f[i - 1], mul(i - 1, f[i - 2]));
		Mul(res, f[c[1]]);

		rep (l, 2, n) if ((l & 1)) {
			f[0] = 1, f[1] = 2;
			rep (i, 2, c[l]) f[i] = pls(mul(2, f[i - 1]), mul(mul(i - 1, f[i - 2]), l));
			Mul(res, f[c[l]]);
		}

		rep (l, 2, n) if (!(l & 1)) {
			f[0] = 1, f[1] = 1;
			rep (i, 2, c[l]) f[i] = pls(f[i - 1], mul(mul(i - 1, f[i - 2]), l));
			Mul(res, f[c[l]]);
		}

		printf("%d\n", mul(ans, res));
	}
}

int main() {

	n = ri();
	rep (i, 1, n)
		G::e.link(i, to[i] = ri());

	G::gao();
	L::solve();

	return 0;
}


Measured memory_bytes: 8704